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PREFACE 


acumen^  is  the  result  of  2  years  of  independent  research  and  study  in 
Cttotf  the  area  of  minimum  distance  theory  and  its  application  to  geodesic  measure¬ 
ment  and  surface  navigation.  The  purpose  is  to  develop  a  new  approach  to 
solving  the  problem  of  getting  from  point  a  to  point  b  on  well  behaved  sur¬ 
faces.  The  study  is  confined  to  the  use  of  curvilinear  coordinates  (u'j1,  uj, 

u3)  and  to  surfaces  which  can  be  represented  by  a  function,  k,  such  that: 

A  .  /-I  r* 

u3  =  k(uA,  i^),  lie,  surfaces  where  one  coordinate  can  be  written  as  a 
function  of  the  remaining  two). 


In  rederiving  the  relations  which  determine  the  geodesics  or  curves  of  minimum 
distance,  we  shall  see  how  using  parametric  separation  of  variables,  along 
with  the  above  coordinate  variables  and  surface  restrictions,  greatly  simpli¬ 
fies  the  problem,  thus  giving  manageable  solutions. 

The  paper  addresses  the  development  of  solutions  not  just  from  the  point  of 
theoretical  interest,  but  in  terms  of  product  useability,  for  eventual  soft¬ 
ware  implementation.  By  product  useability,  I  mean  we  need  to  know  more  than 
just  the  relations  between  u£  and  u3  for  an  arbitrary  minima  curve,  we  need 

explicit  relations  for  distance  traveled  and.  direction  (measured  normally  on 
the  surface)  of  the  minima  curve.  Thus,  we  consider  a  closed  set  of  equations 
the  solution  for  a  given  surface  (c|nd  initial  and/or  boundary  conditions)  to 
be  of  the  form:  l  . 

^  /  7/fP]  ^ - 

u3  =  kCu^,  u2)  (surface  constraint). 


ux  *  f(X)  (parametric  separation,  first  variable). 

u2  3  g(*»  3,  b)  (parametric  minima  solution,  second  variable). 


As  =  s(X,  a,  b)  (parametric  distance  traveled). 

4*  *  'MX,  a,  b)  (azimuthal  heading,  0-360,  off  a  known  constant/known 
direction) . 

where  X  is  the  parametric  separation  variable,  and  a  and  b  are  constants  of 
integrations  which  (assumably)  uniquely  determine  the  curve. 

Several  surfaces  are  addressed,  but  most  rigorously  treated  is  the  sphere, 
where  the  second  order,  nonlinear  differential  equation  which  governs  the 
extremum  solution  is  solved  for  both  boundary  condition  and  initial  condition 
cases,  as  both  cases  are  of  physical  interest. 


THIS  TECHNICAL  NOTE  IS  THE  WORK  OF  OR.  K.  CARLTON-WIPPERN  AND  DOES  NOT 
NECESSARILY  REPRESENT  THE  OFFICIAL  VIEW  OF  THE  US  AIR  FORCE. 
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The  nature  of  the  oblateness  of  the  earth  is  very  briefly  discussed,  and  a 

first  order,  perturbation  approach  (based  on  e  ,  the  square  of  the  eccentric- 

2 

ity)  set  of  equations  are  developed.  A  comparison  of  the  e  perturbed  solu¬ 
tions  to  those  of  the  unperturbed  spherical  solutions  is  discussed  in  appendix 
2 

I.  The  e  perturbed  solution  set  constitutes  the  mathematical  foundation  for 
Project  Frost  (First  order  formulae  relating  to  oblate  spheriod  trigonometry) . 

While  the  purpose  of  this,  document  is  to  provide  useable  products,  not  just 
address  theoretical  development,  some  attention  has  been  paid  to  the  existence 
and  uniqueness  of  solutions,  particularly  in  the  spherical  case,  as  there  do 
exist  traps  for  the  unwary. 

The  application  for  the  equations  and  solutions  for  the  problems  addressed  in 
this  document  runs  the  gamut  from  commercial  to  governmental,  particularly 
mapping  agencies,  and  tactical  and  strategic  defense  organizations. 

The  mathematics  developed  herein  are  for  use  only  as  solutions  to  problems  of 
surface  geometry,  curves  of  minimum  distance.  No  attempt  has  been  made  to 
correct  for  physical  phenomena,  such  as  wind,  geopotential,  tides,  currents 
etc.  Corrections  for  heading  and  positioning  are  outside  the  scope  of  this 
study. 

The  format  used  is  necessarily  terse,  and  the  document  assumes  the  reader  to 
be  proficient  in  variational  calculus,  differential  equations  and  linear 
algebra,  as  well  as  spherical  trigonometry  and  in  other  applied 
mathematics/physics  curricula. 
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I.  INTRODUCTION 

The  investigation  begins  with  the  general  distance  equation: 
(ds,ds)  =  ga^dxadxlD 
or 


for  a  space  curve  given  by  x1  =  x1(t). 

The  curves  of  minimum  distance  (or  geodesics)  are  found  by  solving: 


.2  X  .  a  3 

d  x  +  pX  dx  dx 

I  ds  ds 


®  0 


K&  8g 


3x 


3x 


a£ 

a 


,  called  the  Christoffel  symbols  of  the  second 


kind.1 

For  the  purpose  of  this  investigation,  consider  curvilinear  coordinates  only 

1  2 

(which  greatly  simplifies  things,  reference  Lass  or  Morse  and  Feshbach  ) 

2 

also,  the  investigation  primarily  works  with  R  ,  and  assumes  if  u^  u2  and 
u3  are  curvilinear  coordinates  there  exists  a  surface  function,  k,  such  that: 


u3  *  k(uL,  u2)  .. 


The  distance  formula  becomes. 


^ass.  Vector  and  Tensor  Analysis,  p.  288. 

lLass,  op.cit.,  p.  281. 

2 

Morse  and  Feshbach,  Methods  of  Theoretical  Physics,  p.  47. 


(1) 
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where  hj,  h2,  and  h3  can  be  called  scaling  factors. 

Seek  parametric  solutions  for  Uj,  u2  which  will  minimize  s,  since 
u3  =  ktUj,  . 

The  solutions  (for  Uj.,  u2)  can  be  written  parametrically  of  the  form: 
ux  -  f  (X )  , 


u2  =  g  (X  )  . 

f  can  be  chosen  arbitrarily  (to  facilitate  integration),  find  g(A)  via  the 

4 

usual  methods  of  variational  calculus. 


Let  ii2  -  j  h?(^i) 

1=1 


with  u  =  k ( f (X ) ,  u  ),  u  =  f(A),  then  extremum  is  when: 

0  4  1 


'd  /  35 
dT 


86  _  n  '  _  duZ 

xu :)/  3u*  -2'ljX 


letting  u2  =  g. 


3  dk  8 

T  W  Tg 


(mO'ar  pr  [^(nf) 


x  w2  dk 

+  h3  *  -ar 


8  (  dk  Y 


3g 


=  0 


(2) 


Rearranging  2a,  taking  a  degree  of  freedom  to  facilitate  integration,  for 
example,  assume  there  exists  a  function,  f  (=  Uj),  which  sets 


d 

3T 


l/,2  dk  3  ( dk  \Y]  1  ,2  dk  8  7  dk  ^  n 
rvh3  •  dT  *  8gVdXyJJ  "  I*  h3  *  dT*  8gVdxJ=  0 


(2a) 


3 

Spiegel,  Advanced  Mathematics  for  Engineers  and  Scientists,  p.  288. 
Courant-Hilbert,  Methods  of  Mathematical  Physics,  p.  164. 


2 


1  Ma?ch  1988 


TECHNICAL  NOTE 


then  taking  this  degree  of  freedom  for  f  leaves 


(2b) 


2  -■ 

-  c6 


J*  '  2  2 

n„g  =  c 


but 


(3) 


dk  _  3k_  3k^  - 

dX  =  3X  3g  g 


.4  -2  2.2  -2 

h2g  -  C  h2g 


.2  (  df  V 

2  2  " 

‘  hl  \dX  ) 

V 

+  h3 

Rearranging  and  collecting  terms: 
-2 


Ag  “  +  2Bg  +  C  =  0 
A(X,g)  -  c2(h2  +  h*  (||)  )-  h 


(4) 


,BU,g)=c2h2(§)(|| 


C(X,g)  »  cz 
Now  calling 


(§)(t) 


n(x,g)  -  -  1|b  ±  [b2  -  ac]1/2J. 


And  assuming  n  can  be  uniquely  determined  by  the  boundary,  surface/other  logi¬ 
cal  elimination  conditions,  then 


g  ■  n  (X,g)  .. 


(5) 


The  geodesic  curves  are  found  by  choosing  f  to  make  equation  3  possible,  then 
solving  the  integral  equation: 


gCX) 


'fnlx- 


g)  dX  . 


(6) 


This  can  be  considered  a  navigation  kernal,  in  that  for  a  surface  defined  by 
k,  the  constraint  condition  for  f,  equation  (2a)  and  the  solution  to  equation 
(6)  gives  us  the  parametric  equations  of  the  path(s)  of  minimum  distance,  sub¬ 
ject  to  the  boundary  conditions  (end  points  of  the  curve  on  the  surface). 
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While  traveling  oh  a  curve  of  minimum  distance,  it  is  convenient  to  know  in 
what  direction  to  go  to  remain  oh  the  curve.  'Using  the  natural  definition  of 

the  tangent  vector  •^■^aids  getting  around  on  the  curve.  If, 
ux  =  f  (A) 


u2  *  g(X) 

u3  =  k (X) ,  and  if  &2,  §3  are  curvilinear  unit  vectors,  then: 

*■  3f  3g  3k 

0  =  hl  3X  S1  +  h2  3X  §2  +  h3  3X  §3  * 


(Compare  this  result  to  the  contravariant  vector.5) 


If  the  surface  on  which  the  navigation  is  taking  place  can  be  written  as 
S(uj,  u2,  u, )  =  0,  and  the  parametric  curve  of  minimum  distance  can  be 

written  also  as  a  surface  s(uj,  u2,  u3)  *  0  such  that  the  two  surfaces  inter- 

sect  in  the  space  curve  of  minimum  distance,  then  the  heading  vector  H  can  be 
written  as: 


■4.  -v  -v 

H  =  V?  x  VS  . 

To  demonstrate  this,  let  S,  the  surface  on  which  the  navigation  is  taking 
place,  be  represented  as  locally  (approximately)  planar. 


Similarly,  let  5,  the  surface  associated  with  the  minimum  distance  constraint, 
be  approximately  planar  (this  is  roughly  equivalent  to  using  affine  approxima¬ 
tions  to  the  actual  surfaces  about  a  small,  differential  area  of  local  inter¬ 
section,  assuming  the  intersection  exists). 


5 


Lass, 


op.cit.,  pgs. 


270-272. 
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The  intersection  (of  the  two  planes)  is  a  line.  Let  be  the  surface  normal 
to  S  and  be  the  surface  normal  to  C.  the  direction  of  ft  (if  not  =  0)  is 

clearly  ft  =  ^  x  ft2  by  construction.  In  the  limit  that  the  differential  planar 
surface  approximations  go.  to  zero  ft2+  $S,  ft  -*-$?,  so  that  ft  =  x  ?S 
(gradients  to  be  evaluated  at  a  point  of  interest  on  S) . 

So  long  as  ft  =-  0,  the  magnitude  of  |H|  is  relatively  unimportant.  It  is  the 
direction  of  ft  which  matters.  Eventually,  |ft|  drops  out  of  the  analysis 
entirely,  it  is  from  the  components  of  ft  that  heading  angles  are  determined. 
(The- same  argument  is  true  for  5(A).) 

II  EXAMPLE  IN  CARTESIAN  COORDINATES 

a.  Application  to  the  planar  surface. 

For  this  example  the  surface  is  defined  as  a  plane. 

,z  *  a  +  bx  +  cy 

and  search  for  g(A)  such  that  the  integral  is  extrema.  For  example,  let, 
x  =  A 

y  «  g(A) 

z  »  a  +  bA  +  eg  (A)  , 

then 


5 
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A  + 


(  Z  )2  +  ( b  +  eg"  ) 


1/2 


dX 


=  & 


dX 


36 


With  6  =  integrand,  note  that  =  Oj  therefore, 
— -(  2-~Y»  0,  and  after  integrating  once: 


dX  V  3g- 


r) 


e  =  - 


172 


(e  is  a  constant  of  integration) 


|_1  +  g"2  +  (  b  +  eg"  )  2] ' 

Continuing  with: 

[l  +  g^2  +  (  b  +  cg")2J  e2  «  [g"  +  c  (  b  +  eg"  )  J2  . 

This  is  a  quadratic,  equation  in  g  ,  of  the  form: 

>2  * 

Ag  +  Bg  +  C  *  0  (A,  B  and  C  are  constants) . 

Now,  solving  the  equation  for  g",  let 

g"  *  0,  assuming  0  to  be  uniquely  determinable  from  the  geometry, 
then 


g(X)  -  OX  +  E 

and  substituting  back 

x  ■  X 


y  «  DX  +  E 

z  ■  FX  +  G  [*  a  +  bX  +  c(DX  +  E) )  , 

.  3 

This  is  the  parametric  equation  for  a  line  m  R  , 

If  the  line  passes  through  points  (xQ,  yQ,  zQ  )  and^x^  y  ,  z  ^ 


then, 

x  *  (xi  -  X„)x  *  x0 

y  -  (  Vi  -  y0 )  x  *  y0 
z  •  (zi  '  zo)x  *  ^ 
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A  '  ♦' 

is  an  appropriate  expression. 

Letting  a  =  xx  -  xQ,  b  =  y1  -  yQ,  c  =  zx  -  zQ 

then, 

D  =  al  +  bj  +  cfc* 
and  for  distance 

As (X )  =  (a2  +  b2  +  c2)  1/2  AX  . 

The  parametric  equation  can  also  be  'written 


or  solving  individually  for  X,  then  setting  components  equal  to  each  other 

x  ~  x  y  -  y„  z  ~  z 
0  J  J0  0 
■■■■■■ .  ■  - - ...  M  .... - , 

a  D  c 

The  forms  are  equivalent. 


*0,  a  tangential  vector  (  the  tangential  vector  is  often  referred  to  as 

(jp  -*■  Q  \ 

T  *  so  that  T  *  —  j  J  ,  this  is  equivalent  to  saying  *  constant  on  the 

reference  surface,  where  ^  is  a  heading  angle  measured  on  the  surface, 
referenced  from  seme  known  directional  vector. 
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b.  Minimum  distance  solution. 
Consider  the  following  diagram 


To  find  the  minimum  distance  from  p  to  the  curve  g,  let  £2  be  a  plane  (in  par¬ 
ticular,  the  x,  y  plane) .  If  f  is  the  solution  to  the  minimum  distance 
problan  on  a  plane,  f  is  a  line  of  the  form  y  =  ax  +  b  or  x  =  a'y  +  b^  , 
then  the  minimum  distance  between  p  and  g  lies  along  a  curve  f,  hence  if 
p  =(x0,  y0  )  ,  letting  x  be  the  independent  variable,  then, 

(  is)2  -  (x  -  x0)2  *  (  g(x)  -  y0  )2 

ds 

and  this  is  extremum  when  *  0. 


)  *  (  s(x)  -  Yo) 


dx 


(?) 


so  that 


(sW  -  yt) 


(8) 


Solving  for  x  (on  g)  in  equation  (7),  let  the  result  be  xl#  then  find  v 1  = 

g(  xx)  jassuming  g(  xl  )  *  yQ  and,  constructing  a  line  from(  xQ,  y-Q)  to 

(  xl  ’  yl  )  S0^V8  ^or  As*  Equation  (8)  indicates  that  the  line  from(  xQ,  yQ  )to 

(  xl'  yi  )  intersects  g  at  right  angles.  There  may  be  more  than  one  solution 

to  equation  (7),  in  which  case  the  absolute  (instead  of  relative)  minimum 
would  have  to  be  found,  assuming  one  existed.  Whatever  the  case,  the  solu¬ 
tion,  by  equation  (8),  guarantees  that  the  line  segment  intersects  g  as  a  per¬ 
pendicular. 


TECHNICAL  NOTE 


1  March  1988 


In  that  a  smooth,  well  behaved  surface  can  be  represented  as  a  limiting  number 
of  planes,  (ie,  a  surface  S  can  be  approximated  by  n,  limn+oe  differential 

planar  surfaces  or  affine  surface  approximations)  it  seems  plausable  that  if  g 

A  A 

=  f,  f  a  minima  curve  on  a  surface  S,  then  if  p  is  a  point  on  the  surface  but 

A  A 

not  on  f,  and  if  f  is  another  minima  curve  which  goes  through  p  and  intersects 

A  A  #  A 

f  at  right  angles,  then  f  is  the  curve  of  minimum  distance  between  p  and  f. 

A  <* 

As  will  be  shown  in  the  spherical  development,  f  need  not  necessarily  be 
unique  for  the  general  case.  Stated  mathematically 

Let  p  e  f 

and  let  $  =  £(A) 


and 


*  -  f  (X)* 


Then  $  minimizes  J ds  between  p  and  ■?  if  there  exists  A  such  that 
and 


D($)  •  O(^)  *  0  (f,  f  evaluated  at  AQ  ) 

Now  with  0(f)  =*  h  e  +  7TT  ^  +  h  e 
1  dA  1  2  dA  2  3  dA  3 


hence: 


2  da  da  ,2  db  db  __  ,2  dc  dc 
'l  dA  *  dA  h2  dT  *  HT  h3.  BT  *  eft" 


0 


(apriori  uniqueness  not  guaranteed.) 


a'(A)\ 
b'(A)  . 

c'  (A)  J 


It  must  be  remembered  that  this  case  is  concerned  with  solutions  to  second 
order  differential  equations  with  boundary  conditions.  Oependi  ig  on  the  sur¬ 
face  geometry  and  the  boundary  conditions,  the  minima  solutions  may  not  be 
unique,  as  is  shown  in  the  spherical  section. 


The  Cartesian  example  is  finished  by  working  out  explicity  the  point  to  line 
on  a  plane  example. 

**Compare  this  with  the  definition  of  the  contravariant  vector,  note  also 
ds  ■  (D  •  D)1^2dA. 
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(xq*  yo) 


.y  =  tax  +  b 


(s)2  =  minimum  (x^,  y") 


(As  )2  =(xQ  -  x  )2  +  (y0  -  [mx  +  b]  )  2 

dt"  0  *  "(x0  "  x)+  (“m)(yQ  "  lmx  +  bl) 


hence. 


x0  +  myQ  -  mb 

x  - - - - ■,  call  this  point  x'  . 

1  +  m 


i  •  m(x')*b 


%  ♦  m  yQ  ♦  b 
1  +  m2 


Slope  of  intersecting  line  * 


y  -  yr 


x  -  x„ 


mxo  +  m  yo  +  b 
1  +  m2 

xQ  ♦  myQ  -  mb~ 
1  +  m2 


rnxQ  -  my0  +  b 

.  -  “ 

-  x  -m  x  +  my  -  mb 
0  0  0 


NOTE:  Slope  of  y  ■  m,  slope  of  intersecting  line  »  -  — .  OED, 

III  EXAMPLES  IN  CYLINDRICAL  COORDINATES 

a.  Application  to  the  cylindrical  surface. 


♦(£)*♦  ■'(SOT 


dX  . 


For  the  first  example  define  the  surface  constraint  to  be  that  of  a  cylinder. 


r  ■  a 
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Search  for  8  =  g(X)  such  that  the  integral  is  minima. 
Letting  z  =  aX,  obtain 
s  =  .aJ"<SdX 
where. 


Noting  =  0 


then, 

*  0,  and  integrating  once 

dX  \3g  ) 


+  g'Z  ) 


-1/2 


-2  2  2  -2  C 

g  •  o  ♦  o  g  ,  g  -  6  ,  6  .  j  m 

and 

g(X)  =»  0X  +  9f  ,  9(X)  =•  0Q  +  3X  . 

The  solution  in  parametric  format 
r  •  a 
z  *  aX 
9  *  9q  +  0X 

|also  z  -  |-(0  ”®0)  after  eliminating  xj. 

Now, 

s  •  J"(l  *  82)1^2dX  solves  as 

As  *  a  (l  ♦  B2)  1/2  .  AX  |  . 

*0 

Letting  ax  -  iz  ,  a2  «  ip  ,  s3  -  i0 
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D  =  al  +  rtiQ  =  a£  +  a3iQ  . 

Z  o  Z  0 

-*• 

NOTE:  D  =  Constant,  independent  of  X. 


Similarly,  the  surface  function  is  defined  as  S:  (r  -  a) 
distance  function  is  f:  -z  +  J-^9  -  0Q  )  =  0, 


then, 

H  «  7f  •  7S 


-  -1xz  +  f  (f  )  *e  *  + 


7S  ■  l  hence. 


H 


V  X8 


1  0 


-1  0 


3 


1  *  <> 

3  Xz  l0  * 


H  •  i0 

This  results  in  tan  ip  =>  x - -  3. 

a 

H  •  i 

z 


TECHNICAL  NOTE 


0  and  the  minimum 
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TECHNICAL  NOTE 


NOTE:  The  heading  angle  ij»  is  easily  obtained  when  either  D  or  H  can  be 
written  in  only  two  coordinate  directions,  that  is  to  say,  given  curvilinear 
u  ,  u  ,  u  ,  if  the  surface  can  be  written  such  that  u,  =  constant,  \|)  is  easily 

12  3  K 

A  A 

found  by  dot  product  of  0  or  H  by  the  remaining  unit  vectors  e?,  e*  (i,j  *  k). 

If  not,  rotate  0  or  H  as  necessary  to  produce  a  meaningful  result  (reference 
IVb  for  additional  data  on  the  heading  angle) . 

The  five  equation  solution  for  the  cylinder  is: 

r(A )  =  a 

z(A)  =  aA 

0(A)  =  eQ  ±  (3A 


As  (A 


)  -  a(l  *6 2)1/2X  I  1 


i|i(A)  *  tan  1  8  . 

b.  Application  to  the  spherical  surface.  Another  example  in  cylindrical 
coordinates  is  as  follows: 


2  11/2 


dA 


The  surface  constraint  is  defined  to  be  that  of  a  sphere. 


2  2  I"  .  /  2  1\  l/2l 

+  z  ■  a  ,  [z  ■  ±  ( a  -  r  )  J  . 


2  g 

Choosing  r  =  — - —  ,  which  is  obtained  from  the  following  analysis; 

1  +  A2 


From, 


1/2 


and  since  z  *  f(r)  only 


■\f[Kf)2)(&)2*r2(&y 


dA 


211/2 


dA 
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Then  with  r2  +  z2  =  a2 


dz  p 

dr  z 


TECHNICAL  NOTE 


1  + 


(£j 


2  2 
a  . 


2  2 
a  -  r 


222  H0  > 

Letting  r  =  af  ,r  -  a  f  ,  0  -  g  ,  rrr  =  g 


s  =  a 


.(  1  -  f 


1  f'2 .  fy2V/2  *  . 


) 


36 


Define  the  integrand  as  6,  and  by  inspection,  —  =0,  so,  already  integrating 

3g 


36 

once,  — 7-  ■  c 

3g 


fV 


g  =  c 


'2  ,  2 
g  »  c 


1  -f2 


£2  '1 
+  f  g 


1/2 


>  2  —4 

f  f  ,  „-2'2 


1  -  f2 


+  f  g 


Take  a  degree  of  freedom  by  choosing  f  such  that. 


_-2  ■  -4 

f  v 


1  -  f2 

f'«  f2  ( 1  -  f2  )1/2 


df 


f*(l  -f2)1'2 
csc2udu  *  dX 
-cot  u  *  X 

-  (l.-f 2)1/2 


-  *2  ,  -  X2  +  1  ,  f2 

r  r 


dX  ,  let  f  ■  sin  u,  df  ■  cos  udu 


+  X 


2  * 
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TECHNICAL  NOTE 


2  2  2  2  a 

Now,  since  r  *  af,  r  -  a  f  ,  r  - — ,  then  continuing: 


1  +  X 


1 2  2  2 
r 


/  d£  V  =  —a 
W  (1  + 


(i  +  xa7 

2  /  .  ,  \2 


/.dz  \  /  dz  dr  \  _  a _ 

\dX  /  "  Vdr  dX  )  "  ^  +  x2  )3 

s  =  a  J~6  dX. 


2 

6  = 


X 


♦ Gil 


(  1  +  x2  )3  (t  +  x2  )3  1  +  x2  (l+x2)2  1+X 


—  3  0  ,  =  c  (integrating  once). 

3g  3g 


-2  2  2 
g  =  C 


2  -2  /  ,  .  2  \ 
+Cg  (  1  +  X  J  ,  g 


1  -  C 


g 


(i  -eV)1/2 


g  «  6  t  cos  1  BA  . 


Now,  parametrically: 


(l  +  X2)1/2 


(assuming  r  >  0) 


-1 


9  »  90  ±  cos  8X 


z  •  ± 


aX 


(777rF?  ' 

To  eliminate  X,  start  with 
X  «  8  1  cos  (9”®q)  *  3  1  cos  (  9C 


-  e) 


then  choosing. 
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TECHNICAL  NOTE 


continue  with. 


2  2  -2  2  /-  a  \  2 
r  +  r  3  cos  |  8  -  0Q  J  =  a 

/■ 

rearranging  and  multiplying  by  if 


M  cos2  (8  -  8  )  =  S4 


t  2  2  2 , 

(z  =  a  -  r  ) 


hence, 

r(f*}  cos  (  8  ~  e0)  =  • 

Collecting  factors  and  rearranging  to  Cartesian  coordinates,  obtain 

ax  +  by  +  cz  *  0  which  is  the  equation  of  a  plane  confined  to  contain  the 

origin. 

NOTE:  x  =  y  =  z  -  0  always  a  solution. 

c.  Application  to  the  conical  surface'.  Another  example  in  cylindrical 
coordinates. 


•/[($)*  *(£)2  *■"(£ 


211/2 


Define  the  surface  constraint  to  be  that  of  a  right  cone: 


r  *  az  , 


r  ■  aX 


s  *  aj\<2  *  xY2]1/2dX  (y2  =  1  +  ~  'j 


as  before 


36  n  36 
—  *  0  ,  - -  c  , 


then  rearranging: 

,4  -2  2  2  2,2  -2 

Xg  *  c  Y  +  c  X  g 


E  -  Y 


(f2  -0 


(X  »  fc  ,  dX  »  cdf) 


g  -  8  ±  y  sec  f 
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TECHNICAL  NOTE 


hence : 
r  =  aX 


e.  =  90  +  Y  sec-1  § 


IV  EXAMPLES  IN  SPHERICAL  COORDINATES 

a.  Application  to  the  conical  surface. 

(ds)2  =  (dp ) 2  +  p2  (  d<}> )  2  +  p2  sin2<j>  ( d9  )  2 . 

The  surface  constraint  is  defined  to  be  that  of  a  right  cone, 
<J»  ■  constant,  which  can  be  defined  as 

<f>  *  tan  1  a  . 

Letting 
p  *  aX 
and 

.  2  .  2  a2 

sm  <f»  *  Y  * - r- 

1  +  a4 


continue  with: 

.  -aj(  dX 
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=  0  , 


TECHNICAL  NOTE 


2  2,2/  2.2  2  \  '2 
c  =yx  (yx  - c  ; g 

s'  -  *  fr(A2  -  c2)  -1/2 


gft)  -  %  ±  f  sec'1  2i 


(1)  Conic  distance  measurement. 


As  =  a  I  1  + 


.2.2  2 

d  A  -  c 


r-jh 


V  -c2)1/2 


As.,(x2-^)l/2.f(xV-c2)1/2  . 


(2)  Conic  heading  vector  development. 


s  =  4>  ~  tan  a 


f  *  -  ~  sec  Y  (  9  ~  ®q)  +  P  (choosing  the  positive  direction), 


/XP 

■+■  +  I 

H  *  7f  x  Vs  »  1  0 


p  sin  <f> 


u  tan  u  ( u  *  y(9  -  0Q )) 


ac  .  \  o  1  * 

— r-  — ; —  sec  u  tan  u  l  +  -  xQ  . 

4  ,  .  /  P  0  0 

P  sin  <)>  / 
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TECHNICAL  NOTE 


Now,  looking  on  the  conic  surface: 


cos  i|c  *  j  H  | 


_  |u i  ac  sec  u  tan  u 


tan  ip  = . 


P  sin  (j> 
P  sin  <j> 

’ac  sec  u  tan  u 

Now, 

sin  <J>  =  Y , 

and  from  f: 

ac  sec  y  (  8  -  9  ) 

r  > - 7  ' 

Thus, 


(  sin  *  *  |H |  1 


tan  ip 


tan 


hr ' cot  T  ( 9  •  9o)  • 


Now,  to  show  the  equivalence  between  H  and  0, 
p  =  aX 


<p  *  tan  1  a 


0  =  0  +  J  sec'1  ~ 
o  f  c 


obtain, 

°  *  aiA  +  +  P  sin  ^  (^-)i0 


(p  =  aX 


A 

i 

■  P 


) 


choose 


,  sin  <p  =  y) 
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TECHNICAL  NOTE 


0  =  a 


XP  * 


YX  /y2X2 


ITT  xe 


-  i 


=  a  : 


1 

f  tV.  - 1 1 

\  c2  } 

~T7z  1e 

Now, 

0  •  1  =  a  *  cos  ip  |D  | 


\  1/2 

1  *  sin  ip  (□ | 


hence, 
tan  ip  « 


[(t)2  '  \  1/2  ■(*“2  &(6'-  e0)>l)  ■  cot  y(8  -  80)OED. 


t 

r\ 


Notice  when  the  surface  constraint  can  be  written  u^  -  constant,  the  coor¬ 
dinate  system  can  be  considered  surface  normal.  In  such  a  system,  when  [T  is 
evaluated,  only  two  components  are  present.  These  components  are  along 
orthogonal  topocentric  unit  vectors  on  the  reference  surface.  The  five 
equation  solution  for  the  cone  is; 

p(X)  -  aX 
$  (X )  »  tan  *  a 
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t|/(A)  =  tan 

b.  Application  to  the  spherical  surface.  Define  the  surface  to  be  that 
of  a  sphere: 

P  =  a  (constant)  and  search  for  0  =  g(A)  such  that  the  integral  is  a  minima. 

1/2 

dA 

(1)  Solution  by  parametric  separation  of  variables.  The  differen¬ 
tial  equation  for  spherical  integrating  factor  is  obtained  from, 

45  ■  aJ((*)  * si"2  ♦(<&))  dX  • 

Let, 


s  = 


+  p  sin  <j> 


(£): 


(j)  a  f 
Q  a  g 

anp  obtain: 

As  -  aj"(f"2  +  sin2  fg'2  )  1/2  dA  . 
Noting, 


.  o 

3g  U  ' 

then,  integrating  once: 


36 


sin  fg 


»g  (f'2  •  sin2  fg'2  )  1/2 


2  _2  2  .  2  „  '2  .  k  '2 

c  f  +  c  sm  tg  »  sin  fg 


2  f  2 

C  - - -  +  c 

sin^  f 


'2 


•  2  r 

sin  f 


'2 


■  2  £ 
sm  f 


+  c 


•  2  X 

sin  f 


■  § 
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TECHNICAL  NOTE 


Choose  f  such  that: 


_ 1_ 

sin2  f  X2 

df  _  dX_ 

2  2 
sin  f  X 


2  1 
esc  f  df  =  y 


-cot  f  =  - 


f  =  tan  X 


so  that 


f  =  tan  X 


(X  =  tan  <j>) . 


Then, 


dX  1  ♦  X2  ' 


2  X 

sin  $  „  - -  , 

1  +  X^ 


9  *  g(X) 


s  *  a 


where. 


'f  5  dX 


«  -  If  1  +  V2  +  \z(  1  +  X2 


Noting, 


)■'(*)] 


2-,  1/2 


36  n 

S?"  °  ’ 


<&) 
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TECHNICAL  NOTE 


and  integrating  once  gives 
, 2 


1  +  A- 


g  =  c  <S 


A2  \2  '2  '2 

- r  g  =  c 


c  =  constant  of  integration 


,2 


1  +  A 


\  1  +  A  / 

,  4  '2  '2  '2.2  (  .  ,2\  '2 

Ag  =c  +  C  A  (  1  +  A  )  g 

,4  -2  -*2,  2  '2  '2.  4  '2  '2 

Ag  -cAg  -cAg  =  c 

(  .  '2.  \  .  4  '2  '2.2  '2  '2 

(1-C  JA  g  -cAg  =  C 

,  2  '2  f  f  .  '2  \  ,2  •'21  '2 

X  g  [(1-0  )\  -c  J  -  c 

±  c 


(l  +  A2) 


Setting 


211/2  * 


2  1  -  c 


'2 


8 

and 


q  T 


dq  *  dA 
g(q)  »  ±  J*- 


_^3_ 


q(  q2  -  1  ) 

g(q)  ■  0Q  ±  sec  1  q 
-1  A 


1/2 


g(A)  *  0Q  ±  sec  j  (with  9Q  and  8  constants  of  integration). 


Parametrically,  the  curve  is 
p  *  a 

$  -  tan  1  A 
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0  =  8q  ±  sec"1  j 


TECHNICAL  NOTE 


To  gain  further  insight  as  to  the  nature  of  this  curve,  work  backwards  as 
follows:- 


0  -  0O  =  sec'1  j 
X  =  0  sec  (©  -  0Q  ) 


X 


-c _ 

cos  (9  -  0O) 


a  cos  0  +  b  sin  0 
where, 

a  =  -c0  1  cos  0Q 
b  *  c0  1  sin  0Q 
and  substituting  tan  <J>  for  X 
tan  $  -  a  CQS  Q+5  sin  q 

NOTE:  This  equation  relates  <j>  and  0  explicity,  ie,  $  *  f(8)  or  f  1  (4>)  *  0  . 

Continuing  with: 

a  sin  <}>  cos  0  +  b  sin  .sin  0  *  -c  cos  <J>  , 

Multiplying  this  through  by  the  term  p  and  rearranging  gives: 
alp  sin  <j>  cos  0]  +  b[p  sin  <}>  sin  8]  +  clp  cos  ((>]  3  0 

Now,  recognizing  the  relationship’ between  spherical  and  Cartesian  coordinates: 
x  «  p  sin  $  cos  8 
y  «  P  sin  <j>  sin  0 
z  *  p  cos  <j>  . 

This  constraint  equation,  in  Cartesian  coordinates,  becomes: 

ax  +  by  ♦  cz  ■  0  , 
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TECHNICAL  NOTE 


This  is  a  plane  confined  to  pass  through  the  end  points  (boundary  conditions) 
and  the  sphere  origin.  The  solution  for  set  boundary  conditions  is  unique 
provided  the  end  points  are  noncolinear  with  the  origin. 

NOTE:  x  =  y  =  z  =  0,  always  a  solution  to  the  constraint  plane. 

Since  the  solution  is  confined  to  a  plane  (ie,  linear)  and  the  curve  must  be  a 
GREAT  circle  (ie,  includes  the  system  origin)  the  distance  between  end  points, 
s,  can  be  written: 

s  =  a  A  , 

where  A  is  the  angle  between  the  origin  and  vectors  which  point  to  the  end 
points.  These  vectors  can  be  called: 

pQ  =  a  [sin  <t>Q  cos  8Qi  +  sin  <j>Q  sin  8Qj  +  cos  <f>Q£] 
and, 

'  f  .  „  A  .A  A-l 

Px  ■  a  sin  ^  cos  8^  +  sin  ^  sm  8^  +  cos  (j^k  . 

Note, 

|P0  I  ■  Ifj  I  -  a"  , 
and  from  linear  algebra 


9 


therefore: 

cos  A  =«  sin  <J>Q  sin  4^  cos  8Q  cos  8^  +  sin  <j>0  sin  8Q  simj^  sin  8J 

♦  cos  $0  cos  ^ 

cos  A  *  cos  <j>Q  cos  4>1  +  sin  <{>0  sin  <J>1  jcos  8Q  cos  81  +  sin  8Q  sin  0[J 
cos  A  *  cos  cos  <J>1  +  sin  sin  cos  (8Q  -  8^  . 

This  is  exactly  what  Naperian  spherical  trigonometry  (law  of  cosines)  tells 
the  reader  to  expect.  The  spherical  law  of  sines  can  be  recovered  via 
geometrical  construction. 
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TECHNICAL  NOTE 


If  a  boundary  value  problem  is  solved  in  such  a  way  that  a  minima  curve  goes 
through  two  points,  pQ  and  pi;,  then  the  minima  solution  is  a  plane  of  the 
form 

N  *  x  =  0 
where, 

N  =  (a  b  c)  , 


N  is  the  planar  normal  vector  (plane  of  the  great  circle)  of  arbitrary 
magnitude.  Then  given  pQ_,  p,  N  can  always  be  constructed'  (of  arbitrary 

magnitude)  as  N  =  pQ  x  p1  .  If  N  =  0  when  p  =  ±  pQ  ,  the  solution  is 
nonunique.  In  this  case,  there  are  infinite  solutions  which  can  be  obtained 
by  finding  N  such  that  either 

Px  •  N  *  0  or,  equivalently 

p  •  N  *  0 

(1)  Spherical  detail  on  H,  the  heading  vector.  Define  pQ  and  p^ 

as  points  on  a  given  minima  curve.  On  the  sphere,  prerotate  the  coordinate 
system  such  that  p0  is  in  the  same  plane  as  the  zero  longitude  and  refer  to  pj 

in  terms  of  ,  A9  )  , 

A  review  of  spherical  coordinates. 


N 
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TECHNICAL  NOTE 


With  pQ  at  0  *  O',  and  using  unit  vectors: 

-V  “1  /  A  a  a  \ 

P0  8  a  +  yj  +  zk  ; 


-r-^T  =  sin  «j»  i  +  cos  <|>nk 

iPn  I  °  ° 


and  for  arbitrary  p. 


— 1— ;  =  sin  <j>  cos  A0i  +  sin  <j>.  sin  A0j  +  cos  <p  k  . 
IPi  I  1  1  1 


Now, 


p  p  + 

—2—  x  =  i  (~  sin  <j>  sin  A0  cos  $  )  +  j  (sin  <p  cos  A0  cos  d> 

|P0  I  |PX  I  1  0  1  0 


sin  <j>Q  cos  +  k  (sin  sin  ^  sin  A0) 


also. 


Vf 


|Vf| 


Po 


P0  |P0 


since 


A 

Vf  *  —  (xi  +  yj  +  zk)  ,  f  the  spherical  surface  function, 
a* 


and 


sin  *0i  +  cos  <fr06 


Therefore: 


H  *  N  x  Vf* 
N  *  pQ  x  px 


*By  definition,  also  refer  to  geometrical  interpretation,  since  N  ■  V5,  the 
minima  surface,  (ie,  plane). 
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TECHNICAL  NOTE 


and 


H  = 
r 


^sin  cos  A8  cos  -  cos  <j>1  sin  (j>Q  cos  «j>0  ^ 
sin  4)^  sin  A0 

\cos  ^  sin2  <j>0  -  sin  cos  A0  sin  <j>0  cos  <()0  J 


This  vector  now  needs  to  be  rotated  into  topocentric  (SEZ)**  coordinates. 
This  is  done  by: 


“sEZ 


' cos  <j>0  0  -sin  <j>0 


>sin  <j>0  0  cos  <|>0 


calling, 
c  =«  cos  <j>0 

s  *  sin  <|>0 

x  =>  sin  <j>1  cos  A9 

y  3  sin  <j»1  sin  A0 

z  »  cos  ^ 

/  3  2  3  2  \ 

/  xc  -  zsc  -  zs  +  XS  C  \ 


H 


SEZ 


2  2  2  2 
xsc  -  zs  c  +  zcs  -  xc  s 


/ 

/ c2  (xc  -  zs)  +  s2  (xc  -zs)\ 


:  I 


**SEZ,  a  local  right-handed  coordinate  system  which  points  south,  east  and  up 
(or  Z). 
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N 


c 

By  leaving  the  result  as  a  -fraction,  ie,  tan  t  »  ■=•,  all  quadrants  are  readily 

resolved.  In  an  absolute  sense,  equation  10  could  have  just  as  easily  been 
obtained  via  spherical  trigonometry: 
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(11) 
(12) 

(13) 

continuing: 

sin  c  sin  i|»  =  sin  <j>  sin  48  (14) 

and  from  equation  (12) 

sin  c  cos  *  — . .  fcos  -  cos  <j>A  cos  cl 

sm  <j>0  L  1  0  J 

and  substituting  from  equation  (11) 

sin  c  cos  i|>  *  s'i'n'V "  [cos  ^1  “  cos  (  cos  ^0  cos  ^1  +  s*n  s*n  cos  A9  )] 

*  [cos  *i  0  "  cos2  *o)  ’  008  *o  sin  *i  cos  A8  (sin  *o )] 

a  ii^l7  [sin  *o  (cos  *1  Sin  *0)  “  3in-^o(sin  +1  cos  +0  cos  a03 

sin  c  cos  <|>  ■  cos  ^  sin  <j>Q  -  sin  <j>1  cos  cos  48  .  (15) 

<r 

Dividing  equation  (14)  by  equation  (15)  gives: 
sin  c  sin  ^  sin  'sin  48 

n  ■  .  . S  ih  a  _  —  ■-■■■■ .  ■■■■■.. . — — — ■  ■  —  . 

sin  c  cos  <|>  cos  ^  sin  -  sin  ^  cos  <f>Q  cos  48 

The  geometrical  interpretation  of  H  *  N  x  $f  is  determined  from  the  following 
illustrations: 


cos  c  =  cos  <j>0  cos  4» +  sin  $0  sin  ^  cos  48 

cos  *  cos  <j»0  cos  c  +  sin  sin  c  cos  i|> 

sin  i|>  _  sin  48 

sin  ”  sin  c 
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In  this  illustration  curve  l  represents  the  intersection  of  plane  ft  and  the 
spherical  surface.  Clearly  H  must  reside  in  ft,  therefore  be  normal  to  ft 
(ie,  N  •  H  =  0  always). 


From  this  illustration,  ft  is  the  planar  surface  normal  to  the  point  p.  H 
must  also  reside  in  this  plane.  Headings  (ie,  north,  south,  east  and  west) 

are  locally  defined  in  this  plane.  H  must  also  be  normal  to  M(M  •  H  ■  0 
always ) . 

-v 

NOTE:  N  is  the  surface  normal  for  plane  ft,  M  is  the  surface  normal  for  plane 
ft"). 

H  could  be  defined  as  the  intersection  of  ft  and  ft  or: 

H  »  ft  H  ft' 
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Since  the  magnitude  of  H  is  not  a  factor,  just  the  direction,  then  the  fact 
that  ft  •  ft  =  0  and  ft  •  ft  =  0  can  be  used.  From  vector  analysis  a  convenient 
way  of  finding  a  vector  A  which  is  normal  to  both  £  and  £  is: 

X-txB 

then  define  ft, 
ft  =  ft  xft  . 

Now  ft  is  obtained  from  pQ  and 
ft  -  P0  x  Pj  , 

and  ft,  the  surface  normal  vector  is 

ft  -  ff  , 

where, 

f  =  f(x,  y,  z)  which  defines  the  surface. 

Hence, 

ft  =^pQ  x  Pjjx^f  (for  a  sphere). 

$f  to  be  evaluated  at  a  point  in  question. 
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The  space  curve  directional  vector  5  is  found  as  follows.: 


For  first  quadrant  direct: 


p  =  a 


<J>  =  tan  1  A 

9  =  0Q  -  sec  1  (A  decreasing,  hence  -1  •  dA). 

3a )  -  01  ♦  - i— (-l)i  *  (-1 )- - ■  ,  ■  •  - J — - 

P  (l*x2)  ♦  ( 1  *  x?)1/2  , 

r  U2 


2  u/2(_1)i8  ■ 


Measure  from  true  topographic  north  (ie,  -i^  direction) 

6u) "  (TTI^t"1*)  o  *  *2)1/2(  , ’i 1/2  ie 


is 
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NOTE:  From  this  cot  <j>  =  cos  <j>  tan  (  9Q  -  0)  . 


For  ft  , 

ft  *  x  $S 

Now,  from  the  surface  constraint : 
p  -  a  =  0  , 

and  for  the  minimum  constraint: 

3  sec(9Q  -  0)  -  tan  <j>  »  0  (since  in  the  5  derivation,  a  quadrant  in  which 
0  <  8  was  used,  thus  making  the  argument  positive). 


>  sec2  4  i^  +  p 'Yi'n  sec  A0  tan  40 )  i0  |  ^evaluated  at  p  ,  49)]  . 


Then, 


i„  x, 

P  <j> 


ft  -  0  -«££-  3  sec  A8  tan  A8 

1  P  P  sin  <J> 

\  1  0  0 


_  n4  _  8  sec  A0  tan  A0  «■  sec  $  * 

H  U1p  — p  sin  ? - %  - p  x0  * 
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Now,  in  the  local  system 


- -  -  -cot  9 

H  •  iQ 


3  sec  A8  tan  A9 
P  sin  9 

2  X 

sec  9 
P 


and, 

,  .  8  sec  A9  ,  ,  .. 

cot  9  *  — r - r —  cos  9  tan  A9 

tan  9 

but, 

8  sec  A9  ■  tan  9 

cot  9  *  cos  9  tan  A9  QEQ. 
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Finally,  the  relation  for  could  have  been  derived  directly  from  spherical 
trigonometry. 


\ 

\ 


given  b,  a  ,  find  <j». 
b,  a  implies  a,  c,  d,  ty' 

♦  3  *  "  , 
therefore, 
tan  =*  tan  <|^ 
with 


and 

c  +  b  »  - 
2 


sin 

cos 

-cos 

-cos 


*  -  sm. c. 

sin  s 

d  *  0  *  cos  a  cos  s  +  sin  a  sin  s  cos 

a  cos  s  ■  sin  a  sin  s  cos 

'  _  cos  a  cos  s 
v  "  sin  a  sin  s 


f 
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sin  a  _  tan  c 
cos  a  “  cos  a 

since 

cos  s  =  sin  a  cos  c  . 

Continuing  with 


-tan  ij) 


sin  c 
sin  s 


cos  a  cos  s 
sin  a  sin  s 


sin  c 
cos  s 


tan  a  = 


sin  c 


sin  a  cos  c 


c  +  b 


IT 

2 


I 


tan  c 

tan  i|» 

cot 
but 
a  =  <j> 


=  cot  b 

_  cot  b 
cos  a 

=  cos  a  tan  b 


b  -  A6 


cot  =  cos  <j>  tan  A0 
or 


„  sec  <{» 
tan  A9 

cos  s  ■  cos  a  cos  d  +  sin  a  sin  d  cos  c 


tan  t 
NOTE: 


cos  s  =  sin  a  cos  c  for 


For  distance  measurement: 


As  *  a 


1  +  X' 


1  + 


1  +  X 


-  1 


1/2 


dX 
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=  a  f — •  m~2  ' '  -~T '1/2  [~ ^7  "  1  +  1  +  1 

J  i  +  x  II.-  1 J  Le  J 

r  1  ■]  r 

■■JTT?' 7J7>  H 


As  =  a  (  1  +  B" 


Substituting: 
1  +  X2  *  u 


2XdX  =  du 


(i*x2)(x2-b2)iA2 


XdX  . 


X  =  u  -  1 


a  M  +  B 


2\  1/2 


u(u  -a)1'2 


[«-(s2*l)]  . 


Substituting  back  for  u 


As  *  a 


(l  +  32)1/2  .  2  t  -1  u  -(B2  +  l)  1/2  I  R  ♦  X2  -  B2 

2  ( 02  +  1  ) 1/2  B2  +  1  \  _  B2  +  1 


As  ■  a  tan 


2  .2! 1/2  ,  l 


-1  X  -  B 
2 

1  +  B 


The  five  equation  solution  for  the  sphere: 


p  =  a 


<j>  =  tan  1  X 


0  «  60  ±  sec"1  j 


As(X)  *  a  tan 


2  _  a2  \ 1/2  ,  1 


-1  /  X  -  B 
\  1  +  B2 
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From  the  five  equations  for  the  sphere,  combine  the  equations  for  As  and  9, 
eliminating  the  parametric  variable  X. 

9  =  9q±  sec"1  i 


As  =  a  tan 


-1  /  A2  -  B2  V1/2  I  1 
\  1  *  82  /  A„ 


A  M 

From  the  spherical  geodesic  (great  circle  route)  —  =  A<j>,  the  central  swept 

3 

out  angle,  evaluate  from  the  top  of  the  curve  down  one  side  to  X  =  X  . 


north  pole 


ta„  .  f  LLlil  )  1/2 

a  \  i  +  e  / 


|A9  j  =  |9  -  9  |3  sec 


hence. 


X  3  0  sec  A9 
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Note  that  for  a  given  geodesic  8  is  related  -to  the  minimum  colatitude  value 
that  geodesic  curve  will  attain.  If  =  $,  then  8  =  tan  $.  $  is  related 
to  the  inclination  of  the  geodesic  curve  (remember,  the  curve  is  coplanar)  by 

.  it  i 

i  •  7  -  ♦ 
hence, 

8  =  cot  i  , 

(i  ‘  c°s  1 

and 

tan  -  A ty'  )  =  cos  i  tan  (  ^  -  A0^) 


cot  Ai|>  =  cos  i  cot  A0 


tan  A^»  •  cos  i  ■  tan  A0^  . 

This  is  completely  consistent  with  Napier’s  findings  and  is  sometimes  referred 
to  as  the  spherical  law  of  tangents. 

Although  most  of  the  attention  has  been  northern  hemisphere  oriented,  the 
solution,  by  symmetry,  is  equally  valid  for  southern  hemisphere  trajectories 
(-  *  <  X  <  0).  Singularities  occur  at  X  ■  0  (either  north  or  south  pole)  and 
at  X  -*•  «  at  the  equator.  From  a  computational  point  of  view,  equatorial,  near 
equatorial  and  transequatorial  routes  (X  +  ®)  are  inconvenient  to  calculate. 

Defining  the  variable  x  *  j  and  instead  of  X  «,  let  x  +  0.  For  this  case, 

the  formulas  becomes 

p  «  a 

<j>  (x)  »  cot  1  x  *  y  -  tan  1  x 
9  (x)  =*  0q  ±  cos  1  (8x) 


(x ' 
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For  geodesics  which  cross  the  equator  note  that  tan  i|>  =  0  at  the  equator.  The 
geodesics  are  assumed  to  be  smooth,  continuous  curves  so  that  for  a  given  8  in 
the  northern  hemisphere,  the  southern  hemisphere  continuation  of  the  same 
geodesic  will  have  the  save  value  of  8  in  the  absolute  sense,  thus  both  a 
given  geodesic  curve  and  its  slope  are  continuous  across  the  equatorial 
boundary. 

Other  singularities  of  interest  are: 
lim  8  +  0  . 

This  becomes  a  transpolar  trajectory: 

9  =  9Q  -  ^  on  one  side 

9  ■  0Q  +  ^  on  the  other,  (ie,  0(X)  becomes  a  discontinuous  step  function, 
becoming  undefined  at  X  *  0,  at  either  pole.  Remember  X  >  8  always). 

lim  8  *  ®.  This  becomes  an  equatorial  trajectory:  <j>  i  In  that  X  >  8 

this  becomes  computationally  difficult  to  handle.  The  matrix  solution  removes 
this  singularity. 

lim  X  •>  8  . 

(2 

■  ■* 

X2  -  8 

This  can  be  rewritten  as 

(2 

-- 

X2  +  1 

Thus,  as  X  +  3, 

♦  (X  «  8  )  ■  j  . 

For  a  typical  boundary  value  problem,  given  pQ,  p  ,  find  the  geodesic  curve 

which  contains  pQ,  p^.  The  expression  N  3  pQ  x  p1  is  convenient  to  use.  In 
•  *► 

the  case  of  physical  interest,  where  N  =  0,  pQ  »  -p1  there  are  an  infinite 

number  of  geodesics  which  will  satisfy  the  condition.  (This  is  not  generally 
true  in  the  slightly  oblate  case,  as  will  be  seen  in  IVb (2 ) . 
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(2)  Existence  and  uniqueness  of  solutions.  One  can  use  the  idea 
of  breaking  up  the  surface  into  small,  differential  grids,  and  numerically 
find  a  curve  between  end  points  through  the  patchwork  of  grids  which  will 
minimize  the  overall  distance;  but,  in  general,  there  is  no  guarantee  that 

A  A  <• 

this  curve  is  unique.  Neither  for  the  curve  f,  r.or  f  (as  defined  in  the 
closest  distance  planar  section).  This  is  most  dramatically  seen  on  the 
sphere.  Although  there  are  an  infinite  number  of  cases,  only  two  are  cited. 

(a)  Start  point  north  pole,  end  point  south  pole.  The  minima 

curve  exists,  but  is  not  unique.  Take  any  meridian  to  the  south  pole  (f 
exists,  not  unique). 

(b)  Point  p  at  the  north  pole,  the  object  is  to  get  to  the 

equator.  The  intercept  curve  is  any  meridian  to  the  equator  (f  exists  and 
satisfies  all  conditions,  but  is  not  unique). 

Uniqueness  is  highly  surface  and  end  point  dependent.  Solutions  on  the  sphere 
show  a  parallel  of  sorts  with  the  general  solutions  of  a  second  order  dif¬ 
ferential  equation  as  they  should. 

Given  a  point  and  a  slope,  the  minima  curve  is  uniquely  determined.  Given  end 
points,  the  minima  curve  may  not  be  unique  and  there  may  exist  infinite 
solutions. 

Looking  at  the  spherical  case  in  more  detail,  two  cases  (or  subcases  of  par¬ 
ticular  significant  physical  interest)  are- identified. 

Case  1,  given  pQ  and  <ji*»  determine  the  minima  curve  parameters. 

Case  2,  given  pQ  and  p  ,  determine  the  minima  curve  parameters. 

In  each  case,  assume  pQ  is  at  9  *  0,  (If  not,  use  relative  0  coordinates. 
Remember,  the  surface  is  one  of  revolution  and  symmetry  can  be  used.) 

And  in  each  case  it  is  most  convenient  to  work  with  the  solution  planar  normal 
ft,  which  in  case  2  is  easily  determined  via  ft  3  pQ  x  p^  . 

Now,  consider  the  relation  of  ft  with  ^  : 


a  topocentrically  defined  heading  angle  (0  <  t  <  360) 
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A 

i 


NOTE:  sin  8  «  N  • 

Rotate  the  normal  up  to  the  point  pQ  (at  <j>Q)  this  is  a  y  axis  rotation  through 
an  angle  -  <J>Q  so  that 


/sin  *0 

0 

-cos  <{>0  \ 

1  0  \ 

** 

N  - 

0 

1 

0 

-  COS  ’j' 

^cos  ♦0 

0 

sin  j 

^  sin  V  y 

/ 


-  cos 


N(4>0) 


-  COS  ^ 

\  sin  *0 


sin  \ 


sin  <j»  J 


which  must  then  be  rotated  (-0)  about  the  z  axis  for  proper  inertial  (I,3,K) 
orientation  (0  is  defined  here  as  the  longitudeQof  p  ) . 

This  equation  relates  the  two  cases,  as  given  N,  which  is  equal  to  pQ  x 

and  forcing*  |N  |  *  1  ,  then  i|>  at  <j>Q  is  determined  via  the  above  relation, 

and  given  PQ,  ♦  then  N  is  obtained  (remember  to  keep  track  of  0,  which  was 

assumed  to  be  zero). 
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Now  the  angle  0  is  easily  determined  via  sin  0  =  ft  •  for  either  case 
(reference  illustration)*. 

NOTE:  0Q  the  spherical/sphere  case  parameter  is  obtained  from  5  ^  *  it  x  ft, 

it  5 

from  which  (depending  on  quadrant)  0Q  =  ^  +  tan  x  i 


t  •  i 


By  geometry 


^node  P0*11^  along  the  equator  to  where  the  plane  defined  hy  ft  intersects  the 
sphere  equatorial  plane. 

The  angle  0  is  related  to  the  spherical/sphere  case  parameter  3  via  0  = 
tan  0. 

Given  a  curve  (0n,  0)  and  some  arbitrary  point  p  ,  determine  if  the  closest 

point  on  the  curve  to  point  p  can  be  found. 

0 

From, 

tan  <(>  *  0  sec  (0  -  0Q) 

obtain  the  curve  in  the  cartesian  form 
ax  +  by  +  cz  =  0 

(H  •  x  =  0) 

and  let, 

+  ^  A 

P  *xi+yj+zk  . 
a  a  Jad  a 

Now,  note  if 


ax  +  by  +  cz  *  0 
a  a  a 


then  P  is  on  the  curve.  If  this  is  so,  stop  the  process.  If  not,  continue. 
0 


‘Require  |ft|  =  1. 
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Let  |fl|  =  1  ,  |£  |  =  1* 

a 

construct  %,  S  such  that 
%  •  $  *  0 
* 

a 


So: 


^  «  |flg  •  fl|fl 
$=£  -A  . 

a 

By  geometry  §  points  to  the  solution  point  with 


N 


s  -  -  p  -  flirt  , 

3  3 

then  let 

|8|  =  1 

and  -From  this  8  the  solution  is  extracted. 
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NOfE:  If  P  is  in  the  direction  of  ft,  §  =  0,.  there  is  no  unique  solution. 
Also  if  £  •  =  0,  %  =  0,  $  =  P  ,  the  point  is  on  the  curve  (this  is  an 

3  3 

application  of  Gram-Schmidt  orthogonalization) . 

(3)  Solution  via  matrix  methods. 

(a)  Development.  In  that,  if  ft(t)  is  the  position  vector  of 

the  solution  of  the  spherical  case  (say  the  general  case,  where  P(?)  is  a 
vector  function  of  time  measured  from  the  origin)  it  is  easily  shown  that 
for  a  geodesic: 

Jh  £(t)  •  ft(t)  =  a2,  a  constant  for  any  time.. 

2.  There  exists  a  vector  ft  such  that  f\l  •  F^(t)  =  0,  for 

all  t. 

The  two  results  indicate  a  solution  obtainable  in  terms  of  vector  matricies. 

Consider  the  fallowing  diagram. 

z(0 


Lets  say  ft  is  known  in  the  x,  y,  z  systems  and  |ft|  =  1 
*►  * 

N  •  k  *  cos  8 

-N  x  k  =  0  (which  points  in  the  direction  of  x  ) . 
Now 
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Let  5  have  components  ai  +  bj 


hence 


a  =  tan 


-ft) 


=  tan 


$  •  i 


(all  quadrants  resolved). 


The  known  vector  N  contains  the  two  necessary  components  needed  to  exactly 
specify  the  curve  (of  minimum  distance) . 

Frcm  £(t)  •  £(t)  »  constant  assume  a  solution  in  the  ft  •  x  plane  of  the  form: 


cos  tp 


?(t)  *  sin 


where,  if  the  minimum  distance  trajectory  is  one  of  constant  angular  velocity, 

$  =  wt  +  <5,  <S  an  offset^angle.  $  the  internal  swept  out  angle.  This  vector 
referenced  to  the  x"*,  y  ,  z  system  (a  rotation  of  8  about  the  x  axis)  is: 

/I  0  0  \  /cos  \ 


/I  0  0  \  /cos 

P  (t)  0  cos  8  -sin  8  sin  <J> 


0  sin  8  cos  8/  \0 


and  finally  referenced  to  the  x,  y,  ,z  system. 
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let, 

sin  ♦  cos  £ 
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i 

i 

t 

i 

i 

I 

i 

I 

( 

» 

I 

ft 

t 

l 
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and 

cos  f  >  -  sin  ^ 

for  the  case  of  taking  unit  derivatives. 
Then, 


dft 

dt 


/ 


-cos  a  sin  -  sin  a  cos  3  cos  iji 
cos  a  cos  3  cos  ij>  -  sin  a  sin  <j> 
y  sin  3  cos  if 


\ 


I 


for  ft  x  £ 


\po  P1  p2/ 


i(N  P  -  P  N  )  +  j(PN  -  N  P  )  +  k(N  P  -  PN  ) 
12  12  0  2  02  01  01 


With: 

Nq  *  sin  a  sin  3 


Nj,  *  -cos  a  sin  3 


N2  *  icos  3 

Pc  »  cos  a  cos  <f  -  sin  a  cos  3  sin  if 
Pj  =  sin  a  cos  <f  +  cos  a  cos  3  sin  f 
?1  »  sin  3  sin  <f  . 

Evaluating  each  component 
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(ft  x  p)  •  i  =  -cos  a  sin  8  (sin  8  sin  if)  -  cos  8  (sin  a  cos  if 
+  cos  a  cos  8  sin  if) 

=  -sin  g  cos  a  sin  i|>  -  cos  8  sin  a  cos  i|>  -  cos  a  cos  8  sin  iji 

=  -cos  8  sin  a  cos  if  -  cos  a  sin  if  . 

(N  x  p)  •  j  =  cos  8(cos  a  cos  -  sin  a  cos  8  sin  i|>)  -  sin  a  sin  g  sin  ip 

-  2  £ 

=  cos  8  cos  a  cos  ip  -  sin  a  cos  g  sin  ip  -  sin  a  sin  g  sin  ij> 

=  cos  8  cos  a  cos  if  -  sin  a  sin  if  . 

(ft  x  p)  •  k  =  sin  a  sin  8  (sin  a  cos  if  +  cos  a  cos  8  sin  if) 

+  cos  a  sin  8  (cos  a  cos  if  -  sin  a  cos  8  sin  if) 

.  2 

=  sm  a  sin  g  cos  t|>  +  sin  a  cos  a  sin  8  cos  g  sin 
2 

+  cos  a  sm  8  cas  -  cos  a  sin  a  sm  8  cos  8  sin  ^ 

=  sin  g  cos  ij> 

hence, 

dft  >. 

-jr-  =  N  x  £  (as  expected). 


dft  ± 

dt  *  V?  x*s 


Since  ft  =  $5  and  the  surface  gradient  for  a  sphere  3  what  has  been 
shown  is  that 


dft  ♦ 

»  H  which  is  what  would  be  intuitively  expected. 
Continuing  with. 


ft(t) 


(cos  a 
sin  a 
0 


sin  a  0 

cos  a  0 

0  1 


1  0  0 

0  cos  8  sin  8 

0  sin  6  cos  6 


cos  ij» 
sin  if 
0 


(if  3  u)  (t)  ±  5)  . 
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Finally,  the  term  5,  t  =  0  needs  to  be  evaluated. 

Lets  say  the  term  ^(0)  is  known 
then. 


cos  a  -sin  a  O' 

£(0)  -  (  sin  a  cos  a  0 


/ 1  0  0 
0  cos  3  -sin  3 


0  0  1/a  \°  sin  ^  cos  g/3 

from  which  6  is  easily  solved. 

NOTE:  For  the  argument  6  may  be  =  ai (t )  +  6  or  I  =  oj(t)  -  6. 
The  ±  determination  can  always  be  made  from  the  fact  that, 

$ 


cos  6  = 


i  h  1 1*  i 


.  ,  ft  X  5  I 
ft  1 1*  I 


Direct  trajectory  can  be  defined  for  oi  and  retrograde  trajectory  for  -w  as 
usual. 

Now,  look  for  heading  information  for  £(t). 

Given  £(t)  and  with  known  factors  {$,  9}  or  {latitude,  longitude}  data 

—  can  be  found.  Then  with  the  {latitude,  longitude}  data  rotate  this 

vector  ^  into  the  topocentric  SEZ  system  (the  Z  component  of  ~  must  go  to 
zero  in  this  system)  and  obtain: 
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d P'  * 

W  ’  S  =  cos  + 


dr 

dt 


•  E  =  sin  t 


hence  (where  S,  E,  Z  +  V ,  Y,  k) 


dP  £ 

tan  * .  IT  e . 


^  *  -dp'  .  r 

dt  ‘ 

£  is  the  position  vector  in  the  topocentric  system,  related  P  via 
sin  L  0  -cos  L  \  /  cos  8  sin  8  0 

P'(t)  -|  0  1  0  -sin  8  cos  8  0  )  ?(t)  . 

cos  L  0  sin  L  /  \  0  0  0 

cos  L  cos  8 

NOTE:  L  and  8  are  obtained  from  ?(t)  as  £(t)  ■  (  cos  L  sin  6  1, 


sin  L 


Let, 


£  *  sin  $  cos  8i  +  sin  <j>  sin  8j  +  cos 


and,  noting  <j»  +  L  »  — 

2 


(L,  the  latitude  variable) 
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so  that 


^  =  [  -sin  9 


cos  (J)  cos  8  cos  <j>  sin  0  -sin  $ 
cos  9  0 

sin  <j>  cos  9  sin  <j>  sin  0  cos  <j> 


^  =  j  0  I  =  k 
1 

This  is  exactly  as  it  should  be  since  ?  =  k  as  seen  in  the  local  topocentric 
system.  Letting  M  be  the  topocentric  rotation  matrix, 

h  - "  f „ 


Now, 


K 


/ 


sin  <j>  cos  0 


\ 


/  cos  a 


cos  ^  -  sin  a  cos  8  sin 


sin  <)>  sin  9 


*\ 


sin  a  cos  ijJ  +  cos  a  cos  8  sin  if 


cos  <j> 


yM>*9  \ 


sin  8  sin  $ 


/a,  8,^ 


(this  relates  $,  8  to  a,  0,  iji.  ip  is  the  internal  swept  out  angle). 

•  i  ■)> 

For  the  heading  vector  development  at  the  point  P,  let  H  »  H  H.  H  is  then 
the  heading  vector  in  the  topocentric  system.  '  1 

dfc 

HraMar  » 

NOTE:  hi  •  k  h  0,  since  at  the  surface  point  ?i  is  tangential  to  the  surface. 

ty,  the  heading  angle  as  measured  in  the  local  topocentric  plane  can  be  obtained 
as  usual. 


*  tan 


-1 


V  J 

+ 

H  •  i 
r 


(note  the  negative  sign.) 
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(i) 


(b)  Summary  of  matrix  solution  for  the  spherical  case. 

Let  £  be  the  position  vector: 

(sin  <|>  cos  8  \ 
sin  <J>  sin  0 
cos  <)>  / 

From  the  surface  constraint 
1^1  -  a  . 

From  the  minimum  distance  constraint  obtain 

P  •  ft  -  0  . 

cos  £ 

Let  [  sin  £  ]  be  the  vector  measured  in  the  ft  •  %  «  0  plane  (central  angle, 
0 


1 


V,  measured  from  the  equatorial  node).  Then  9Q,  g  initial/boundary  parameters 
can  be  translated  to  a,  if  angle  constraints  such  that 


cos  a  -sin  a  cos  y 

P  »  (  sin  a  cos  a  cos  y 

0  sin  y 


8  ,  g  determine  ft  and, 
0 


sin  a  sin  y  ^ 
-cos  a  sin  y 


cos  y 


/ 
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=  sin  a  sin  Yi  “  cos  a  sin  Yj  +  cos  Yk  . 


^  can  be  written: 


i|/(t)  =  i|>  +  —  At  (for  constant  velocity  trajectory). 

U  3 


cos  a  -sin  a  cos  Y  sin  a  sin  y 

ft  =  (  sin  a  cos  a  cos  Y  -cos  a  sin  y 

0  sin  y  cos  y 


^  (ft,  the  directional  vector) . 


can  be  determined  as 


cos  4/ 


Z  _  Po„  *  (  %  x  ^  \  if  If5  I  -  -] 

°  xfi|  ' 


NOTE:  All  development  was  based  on  the  solution  to  second  order  differential 
equations. 

The  sphere  was  treated  with  both  condition  types,  ie. 

Case  1 : 

OE :  Sphere 

BC:  ft(t  *  0)  =«  ftQ  ,  ft(t  ■  t)  =*  £  (boundary  value  problem). 


Case  2: 

DE :  Sphere 


dft 


IC:  ft(t  *  0)  ■  ftQ  »  dF  ^  *  0)  *  tan  ♦  (initial  condition  problem). 

In  both  cases,  it  was  convenient  to  find  ft,  the  planar  normal. 

Both  cases  have  physical  interest  since  case  1  is  the  classical  navigation 
problem,  case  2  can  be  thought  as:  given  the  position  and  direction  of  an 
enemy  aircraft,  where  is  it  going,  what  will  it  most  likely  attack,  minima 
curve  assumed. 

In  the  second  case,  an  uncertainty  in  (ie,  surface  directional  heading, 
where  ^  ±  AE)  will  generate  a  solution  envelope. 

The  solution  envelope  is  easily  generated  if  ft.  is  rotated  to  the  north  pole. 
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The  solution  obtained  is  a  surface  envelope  (or  sectoral  slice)  which  can  be 
useful  in  determining  probable  targets. 

(4)  Selected  problems. 

(a)  Intersection  of  minima  distance  curves.  In  the  plane,  the 
following  expressions  can  be  written  for  the  two  minima  trajectories: 

y  *  ax  +  b 

=•  cx^  +  d 

The  objective  is  to  find  the  minima  trajectories  when 


*1 

1 

ii 

i 

i! 


y  *  /  , 

X  *  x' 
hence, 

ax  +  b  *  cx  +  d 
(a  -  c)x  *  d  -  b 
therefore. 


x  *  x  * 


d  -  b 
a  -  c 


(a  *  c) 


y'y'  ■ a  (fMf)  * b  ■ c  (fHr)  * d 

'  -(rrs) 

(rh) 


(ad  -  ab  +  ab  -  cd)  = 


a  -  c 


(ad  -  cb) 


(cd  -  cb  +  da  -  cd)  ■ 


a  -  c 


(ad  -  cb)  (results  are  consistent). 
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NOTE:  If  a  =  c  there  is:  no  solution  if  b  *  dj  no  unique  solution  if  b  =  d. 

Now,  no  unique  solution  implies  a  *  c  and  b  =  d,  ie,  the  trajectories  are  on 
the  same  curve. 

A 

No  solution  implies  y  and  y  are  nonintersecting  parallel  lines. 

On  the  sphere  there  exists  similar  cases. 

Let 


ft  •  x  =  0  (one  curve) . 
fiT  •  x  *  0  (another  curve) . 

The  vector  which  points  to  the  intersection  of  the  curves  is  easily  seen  to 
be: 

NOTE:  The  intersection  vector  t  must  reside  in  the  plane  ft,  hence: 

1  •  ft  *  0  ,  but  1  must  also  reside  in  ft”*  so 

1  •  ft  a  0  ,  so  that  I  is  orthogonal  to  both  ft,  ft"  therefore 
1  *  ft  x  ft  ,  now  since  |f |  »  a  (if  1  *  0),  t  can  be  set  equal  to  ft  x  ft  ,  then 
let,  |t|  -  1.  t  now  points  to  the  latitude/ longitude  of  the  intersection/ 
intersections  (remember  there  are  two,  1  and  -t) . 

£  *  ft  x  ft"  (since  I  must  be  a  vector  in  both  planes). 

1  is  unique  if  and  only  if 

ft  x  ft"  ^  0  . 


If  ft  x  ft  ■  0,  then  either  ft  ■  ft"  or  ft  »  -ft"  which  is  similar  to  case  2  in 
the  plane,  no  unique  solution.* 

There  is  no  analog  to  planar  case  1  (no  solution)  on  the  sphere. 

(b)  Trajectory/intercept.  In  Cartesian  R  let  the  invader 

parameters  be. 


V*> 


♦V 


9 


and  the  interceptor  parameters  be 


*The  trajectories  are  on  the  same  planar/sphere  intersection  curve. 
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^(t)  =  ^  +  t/xt 

where  %  are  assumed  constants. 

Then  at  some  time  (t  =  x)  it  is  required  that 

fyi)  =^(T)  . 


The  object  is  to  solve  for  the  components  of  given  ^  and  . 
Let, 


\  -  V  V  z0 
"  V  V  Z1 


v  =  v  ,  v  ,  v 

0  xO'  yO'  20 


V,  ■  V  .  ,  V  .  ,  V  . 
1  xl  yl  2l 


(Ki  *^i  +  Vzll  -  and  it  is  assumed  |  known^ . 


Now  at  t  =  x 


*0  *  VX0T  *  *1  *  VxlT 


yo  *  V  •  yi  *  V 

*0  *  V  •  21  *  VzlT 

isolating  the  V1  terms: 
Vxi  ■  7  (  *o  •  *1  *  VxoT 
Vvi  ”  r  (  yo  ‘  yi  *  \ 


V  ■  -  ( 
zl  x  V 


yo 


*0  '  21  *  VzlT 


) 

) 

) 


(16) 
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*•  »/ 

hence  we  have  a  quadratic  equation  in  x 

'fr*  -  K  '  *1  *  VxOT  )2  *  (y0  '  yl  *  vyOT  ^  *  (20  "  21  *  VzOT  )2 
^2  -(><0  -  xi)2  *(yo  -  yi  )2  *(zo  -  zi)2  *  vot2  *  2(x0  -  Xl)UxOT 

*  2  (yo  •  yl  )  V  *  2  (20  -  2i  )VZ0T 

K  -4  y1  -  |p0  -?ii2  *  2(^o  -  M*  v 

f  2  2  \ 

or,  assuming  1  >  VQ  )  this  can  be  rewritten  as : 

2 

at  -  2br  -  c  =  0 

b-(fo  M'?o 

-■If.-M2 

which  has  solutions: 
x  =  1  [b  ±(b2  +  ac  )1/2]  . 

Assume  x  >  0,  so  the  solution  is  physically  real,  therefore: 

1  -K  -iVkh  -  D(po  -pi  )•  5»i2 

*  -\i2K  -vo)]1/2] 

The  components  for  V ^  are  found  via  substituting  x  into  equation  (16). 
Thus,  the  values  of  both  ^  and  x,  the  time  of  intercept  are  obtained. 


Now 


if  ( <  vjj  ]  ,  a  physically  real  solution  can  only  exist 


if(fo  -fi)  •  «.  <  0 


For  example,  let 
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then  , 


t  * 


P.  -V 


-2l\  ) 


t 


T  T  2 

(since  |P0  -  P1  >0,  then  for  a  solution  to  exist 


(VM- V° 


(The  case  of  V 


2 

0 


the  discriminant  must  he  real  for  a  solution  to  exist.) 


(c)  Curve  of  pursuit. 

NOTE:  Although  this  section  works  the  curve  of  pursuit  problem,  the  intercep¬ 
tion  problem  on  a  sphere  is  also  addressed  here.  This  is  because  the  solution 
technique  is  the  same  (or  very  similar)  for  both  problems. 


.•(*!»  yx) 


Now, 


ie,  let 

h  ■  *10  *  V 


*1  ■  *10  *  V 


and  let 


VQ  cos  8 


Voy  *  V0  Sin  9 

then  set  a  parametric  relation  for  x(t),  y(t). 
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#•  «  W 

$(t)  •  v 

x  (t)  -  V 

r 

where: 

2  2  2 
r4  =  (x  -  x)  +  (y  -  y) 

V  =  |V|,  the  velocity  (assumed  constant)  of  the  interceptor. 

The  idea  is  the  instantaneous  slope  (or  direction)  of  the  interceptor  points 
along  the  line  (curve  of  minimum  distance)  which  includes  both  interceptor  and 
target  instantaneous  positions.  As  with  most  differential  equations  of  this 

type  exact  solutions  depend  on  (x(t),  y(t)l  the  target  instantaneous  position 
and  the  initial  conditions. 


On  the  sphere,  ft  refers  to  the  target  and  ft  to  the  interceptor.  As  ft  of  the 
target  changes,  iterate  on  ft  and  ft  of  the  interceptor.  It  is  convenient 

(for  ft  =  constant  target  on  a  minima  course)  to  use  the  matrix  concepts 
developed  earlier.  Treat  this  case  analgous  to  target 

y  =  y0  +  vyt 

x  >  xQ  +  vt  . 


The  instantaneous  ft  is  always  ft  x  ft.  This  is  used  as  the  motion  in  At  is 
divided  up. 


Say  for  example  the  given  terms  are:  ftfl,  ft,  (|V|, |V|),  ft,  ftQ  *  x  ^0  '  ^ien 


calculate: 


/ 


P1  »  ft 


\ 


□s  ^6  +  1^1  At^ 
in  ^  At^ 


\ 


/ 


-  .  ± 


(remember  M  =  M(N);  M,  related  to  the  a,  3  matrix. 


reference  IVb ( 3 ) ( a ) 2 ] . 


63 


1  March  1988 

and 


TECHNICAL  NOTE 


=  M. 


j  COS  ^6  ±  ^  At  \ 

in  (&  ±  I-l  )  At 

/ 


\ 


M  =  M(fi) 

then  recalculate  ft. 


ft  =  ft  x  ft 
'l  rl  *  rl 


Since  ft  is  constant,  M  is  constant  and 
j  cos  +  n  ^  At^  ^ 

+  n  ^  At^ 


£  »  M 
n 


sin 


\ 


/ 


Now,  since  ft  changes  for  each  At  cycle,  M(ft)  changes  and  for  a  given  cycle 
(call  it  M) 


ft  -ft  xft 
m  m  m 


-  M(ftj 


m 


m 


pm+1  * 


fcos(S m  *  ^  4t)^ 

sin  fs  ±1—1  At^ 

*  \  m  a  J  I 

\  0  / 


and  continue  the  cycle. 

As  in  the  planar  case,  the  spherical  case  is  easily  generalized  to  PQ 
for  P  a  general  curve. 


P 


(nAt) 
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NOTE:  The  £  terms  are  ?  (s  ±  *  At  )  .  The  ±  can  be  resolved  by  requiring 

i^l*  Pj  >  or  ,  so  the  interceptor  always  points 

towards  the  target. 

(Assuming  |V|  >  |V | ) . 

(d)  The  spherical  interceptor  problem. 

Given  (and  \IQ)  and  ^  (and  find  ft  ,  the  intercept  path. 

Let  be  the  point  of  interception. 

For  the  two  vectors  A  and  5  and  a  unitary  matrix  M 

A  •  $  *  mA  •  (since  only  the  relative  orientation  is  changed) . 

Now  let 


cos  ij>. 


£t  =  M  sin  >J>1 


0 

cos  6 


=•  M  (  sin  6 


( "■*>, )  •  R «  o] 


[-• 


the  a,  0  matrix  defined  so  as  to  make 


^  Ai  not  proportional  to  x  ^  has  been  implicity  assumed.  The 

^  X  £  M 

assumption  should  be  apriori  checked  by  — *  —  , 

i*,*M  i*i. 
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Hence 


cos  i|>. 


M  =  |  sin  i|/j 
0 


cos  <5 

/ 

M  1 PQ  =  (  sin  S 


(6  known). 


Continue  with 

I  cos  (6  +  wx)  \ 


M_1I^0(t) 


sin  (6  +  cot) 

\  0 


(-■?)  ■ 


fflfc  =» 


b 

W  / 


(fixed)*. 


Then  find  x  such  that: 


V. 


M  lP0(T)  •  Fi  ^Pg  CO)  *  cos  -p  x  (target  equation) 


and 


V, 


-It  -  l  1 

M  Pj  •  M  P0(x)  *  cos  —  x  (interceptor  equation). 

j  cos  <S 

sin  6  |  ,  so  the  first  equation  is  always  true. 

\  0 

The  second  equation  becomes 


M  1l^0  (0)  is  defined  as 


*PQ  the  target  vector  P^,  the  interceptor  vector. 
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a  cos 


(6  +iLt)  +  b'  sin  (6  +^-T)  - cos  (^0* 


Now  since  (for  physical  solutions)  t  >  0,  the  x  =  f (x)  type  iteration  solution 
can  be  used: 


=  y-  cos  1  Ja  cos  ^6  +  -gL  x^  +  b  sin  ^6  +  ^  T) 


starting  with  x  =  0  and  continuing  until  some  convergence  criteria  is  met. 
Then  with  x  solved, 

(x)  =  P1  . 

With  the  point  of  interception  found  known  )the  interceptor  can  be 
navigated  to  the  target  on  the  ft  =  x  plane. 

(5)  Directional  range  and  location  problem.  ' 

Given  ,  s,  t  find  £  :  (eg,  let  PQ  be  a  radio  beacon  at  a  known  location. 
Given  a  radio  receiver  at  ^  which  can  determine  it’s  distance,  s,  from 
and  the  direction  (i|0  to  PQ  ,  where  is  the  receiver  located,  determine  £,.) 


r 


has  coordinates  <j>Q,  9Q,  and  ^  has  coordinates  $  ,  0  j  $0,  0Q,  >J>  and  s  are 
known. 
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Now 

sin  i|<  =  sin  A9 

sin  <j>Q  sin  s 


so  that 

0l  =  0O  +  sin  1  (sin  s  (0l  determined). 

NOTE:  Quadrant  relationship  to  A0  positive. 

Also: 

cos  s  =  cos  <j>0  cos  ^  +  sin  <j>Q  sin  <f> ^  cos  A0 


cos  <j>0  =  cos  <j>1  cos  s  +  sin  4>1  sin  s  cos 


sin  <j>,  = 


cos  s  -  cos  4>0  cos 


1  sin  <f>0  cos  A0 


cos  =  cos 


s  cos  <f>  +  (  4- v-"— ^  —f-fl  \  (  cos  s  -  cos  <j>  cos  i) 

l  sin  <}>0  cos  A0  J  \  o  T1  / 


sin  s  cos  f 


)/  w>ill  3  y  VJC  V 

-  COS  4> .  COS  S - : - : - — — 

l  \  sin  <J>0  cos  A0 


sin  s  cos  cos  o 


cos  (j). 


cos  <}>  = 


cos  <J>0  sin  <J>0  cos  A0  -  cos  s  sin  s  cos  t 

- : — r - r-z - - - : - -  (  <j>,  determined  ) 

cos  s  sin  <j>g  cos  A0  -  cos  <t>Q  sin  s  cos  i|»  \  1  / 

sin  2<j>  cos  A0  -  sin  2s  cos 


1  2  cos  s  sin  <J>0  cos  A0  -  cot  <j>0  tan  s  cos  ^ 


NOTE:  On  an  oblate  planet,  if  it  can  be  assumed  the  locator  beacon  travels 
the  minimum  distance  route,  finding  P^,  given  P  ,  \p,  s  is  an  iterative 
process.  Start  with  the  sphere  relations  and  iterate  using  oblate  geodesic 
curve  relations.  Solution  is  unique  as  long  as  P  is  not  a  pole  (north  or 
south) . 

cos  s  sin  <J>Q  cos  A0  t  cos  $  sin  s  cos  ^ 


but  if  4>0  *  0  then  cos  »  cos  s,  A0  undeterminable. 


(6)  Given  range,  calculate  location  problem 
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(a)  On  a  plane. 

Given  ^ ,  P2,  jsj,  |s2|  locate  Pfl. 


and  assume  ^  unique,  so  that 

/  /  N  2  2 

( xi  "  xo)  +  (h  -  y0)  =  si 

2  2  2 

(  x2  -  xQ)  +  (  ~  ^0  )  “  s2  *  ^W0  ecluations  in  two  unknowns  (note, 

there  may  be  two  solutions,  one  solution  or  no  solutions). 

Now 

*0  ■  *2 1  [4  ■  -  y«)  T/2 

so  that 

v*2  m±  K  '(ya '  yo)T/2  4  R  - (vi  - y») 2]1/2 

and  solve  for  yQ  (then  for  xQ  )  . 

(b)  On  a  sphere. 
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Let  . 


£  •  ?  =  known  =  d 

0  1  1 

£  «  £  =  known  =  d 

0  2  2 

|?Q  |  -  known  •  1 
hence,  assuming  |?.  |  =  1 


(7)  Triangulation  problem.  Given  and  <{»  ■find  t,  the 

point  of  interception. 
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aQx  +  bQ  =  a]Lx  +  bj 
and 

bi  -  bQ 

(  Provided  a„  *  aj ) 

y  =  aQx  +  bQ  =  ajX  +  b1  (either  relation). 
Hence 


t  solved  for: 


(b)  On  a  sphere. 
Also  trivial,  look: 


TECHNICAL  NOTE 


As  previously  done  (reference  iVb ( 2 ) ( b ) 1  P. ,  ip.  generates  ft.  and  t  is  gener¬ 
ated  as:  11  1 

then  setting  |I|  ■  1  obtain  (<j>,  0}  directly.  Always  make  sure  the  quadrant 
is  correct  as  ftQ  x  points  to  ±  f, 

c.  Application  to  the  slightly  oblate  spheroid.  First  order  perturbation 
techniques  are  used  throughout. 

(1)  Surface  modeling.  The  zero  order  approximation  is  that  of  a 

sphere. 

2  2  2  2 
x  ♦  y  +  z  *  a  or  p  «  a  . 
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Working  with  the  cross  section:  x  +  y  =  a  . 

For  the  first  order  planetary  shape  approximation,  the  assumption  is  made  that 
the  surface  normal  points  roughly  in  the  direction  of  the  acceleration  vector. 

(a)  Nonrotating  surface. 

flj  ,  v(x2  +  y2  -  a2  )  , 
let 


Also, 


R  a  scaler  constant. 


tan  4> 


dx 

dy 


or,  the  y  acceleration  component  divided  by  the  x  component. 

(b)  For  a  rotating  surface.  The  x  component  is  reduced  by  an 
amount  -w  x  due  to  centrifugal  acceleration.  If  the  resulting  surface  is 
still  to  remain  approximately  normal  to  the  acceleration  vector,  then 
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£r 

dx  _  _ R  ^ 

’^'7^  2\ 

V  R  “  /  X 

yielding 

(  «  2)  2  2  ,2/2  oj2R 

V  1  -  e  )  X  +  y  -  b  I  e  =  g 

which  is  of  the  form  of  an  oblate  spheroid. 

Now,  for  a  solid  body 

O'  <  g  <  gR» 

since  g  =  g(r)  (assume  symmetric,  isotropic  body). 

4tt  g  —  — 

g(r)  =  pr  p  =  average  density  (assume  .symmetric,  isotropic  body). 

As  a  first  guess  (and  trying  to  strike  a  compromise  between  a  solid  body  and  a 

4 

thin  surface)  let  g(r)  «  gR*  the  arithmetic  average  of  g(r)  on  r,  from 
0  <  r  <  R,  R  the  approximate  radius  of  the  spheroid  volume. 

Then  substituting  back: 


2irGp 


For  the  earth 

w2  =  5.3176*10  9  (radians/second) 2 

G  *  6.673*  10“H  N  M2/kg2 

P  =  5.522  gem  3 
yielding 

<  e2  >  3  0.00689  . 

The  plumb  line  is  deflected  as  follows: 

1 

*The  same  value  is  obtained  via  <  gp  >  *  R 


R 

f  g(r)  dr. 
0 
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and 

sin  e  _  sin  <j> 

2  v 

«i  R  cos  <j»  &R 


hence 


0) 

2 


sin  2*  B.  -  a  a  a 


2gc 


2  n 

For  an  alternative  derivation  of  <  e  >  ,  dimensional  analysis  can  be  used. 
From  the  Buckingham  pi  theorum  choose  the  following  pi: 

*  e,  eccentricity  (a  pure  number) . 

7r1  «  w,  angular/ rotation  velocity. 

ir2  ■  p,  density. 

ir3  »  G,  gravitational  constant  (although  a  constant,  watch!!). 


Start  as  follows: 


SFowles,  Analytical  Mechanics,  P.  130. 

7 

Hughes,  Schaum's  Outline,  Fluid  Mechanics,  Chapter  5. 
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Assume  ui  in  the  numerator  (this  seems  feasable  as  cot,  et)  then  going  back  to 
the  pi  theorum: 


NOTE:  This  also  assumes  n  >  0. 


-i 

If  n  =  •sf  is  chosen,  then  the  value  for  k  is 


X  /  2 

5  0.89099  (  <  et  >  -  0.08301  )  . 


theoretical 


actual 


s  0.68253  (  <  e  >  s  0.081992  )  . 


/  kt  "  ka  \ 

Thus  giving  a  percent  .ifferential  for  k,  (100  •  - ^ - 1  “  1.22^!! 

(2)  Oblate  spheroid.  According  to  1975  CRC8  for  the  earth, 
a  »  6378.388  km 


b  -  6356.912  km 

yielding 

e2  =  0.00672  . 

NOTE:  These  sections  are  not  meant  to  serve  as  a  rigorous  derivation  for <&, 
the  planetary  eccentricity  value,  but  more  as  a  heuristic  demonstration  on*  how 
to  arrive  at  an  approximate  value. 

Consider  the  surface  model  to  be  a  slightly  oblate  spheroid  with  eccentricity 

2  2 
e  <  1,  e  <<  1.  In  the  analysis  consider  only  terms  up  to  and  including  e  . 

All  other  terms  shall  be  discarded  as  being  vanishingly  small.  Continue  with: 
2  2 

~  +  ~  *  1  (letting  r  *  p  sin  $  and  z  *  p  cos  <j>) . 
a1 


a 

Consolidated  Rubber  Company  (CRC  Handbook  of  Mathematical  Tables*  1975). 
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2  2  2  /  2  .  2  x  2  2  \ 
P  -  a  b  ^  b  sin  <j>  +  a  cos  <j>  ) 

2  .2 1  b2  .  2  A  2 

P  *  b  I  —  sm  <}»  +  cos  $ 


-l 


)" 


V.b1  h-i-LW  ♦)- 


,  .  ..  b  2 
let  1 - -  =  e 


-  b(l  - 


2  .2  n-1/2 

e  sin  <(>  J 


Now,  assume  e  <<  1. 
2 

r* 


0 


P  =  b  ^  1  +  tj—  sin  <j>^ 


and 


d£ 

d<j) 


"  ~  order  of  e 


.4 


set 


dp 

dX 


2  '  0  (sln“  ar  -  if  ax)  • 


Continue  with 


-fb0  *f-sir|Z  ♦)  (ir)  *sin2i»(ir) 


n  i  / 


1/2 


dX 


As  in  the  case  of  the  sphere  let 


<f>  *  tan  1  X 


d*  1 

m »'  l  X  -  r - 


dX  1.x2 


.  2  X‘ 

sin  (|)  » 


1  +  X' 


so 


s  *  b  1  + 


e2X2 


2(1  ■+  X2)  /  \  (1.x2) 


X2  p"2  \  1/2  ri\  f 
—  g  )  ^  ( 


since 


-  d9  \ 

=  HI ) 


77 


1  March  1908 


TECHNICAL  NOTE 


so  that: 


"i ,  i>L  ~ 
.  2(1*  \2)_ 


1  +  X 


r  ,2 

- - - +  ,  -  A  g'2- 

(1  ♦  x2  )2  1  +  x2 


letting 


c  =  c  I  1  + 


e2X2  \"1 


V  2(1  +  X2)  j  • 

c'2ac2/i 

l  (1  ♦  *  )/ 

^  a2  "  \  2  "'2  /  1  X2 

l  -1  +  x2  j  =C  l(i  +  x2)2  \  +x2 


X  <2 

+  -  O' 

9  « 


,4  '2  ^2  'Zy2  ( A  .  2  \  -2 

Xg  ■  c  +cX(1+.\jg 

,4  -2  2  2  2  X2  2,2  /,  ,  2  X  '2  2  2,4  'J 

X  g  -  c  -  c  e  - “  +  oX  \1  +  A  jg  -ceXg 

1  +  X2 

(1  *  cVlx'g'2  ■  c2f-1  "  j-  c2X2(1  *  X2)g'2 
x2rx2n  *  cVi  -  =2(i  *  x2)V2 .  c2/ 1  -  -&l  \ 
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f4-iV2-(i-^y 

?  J  \  1  +  A  / 


1  1  /  .  2  2  2  \ 
-r  a  ~7  \1  -  0  +  e  c  ) 

r  c 


[\2  11/2  - 
V  - 1  g  =c± 
F 


/  2y2  \  2 

)  1 - - 2\  (again  expanding  on  e  ) . 

\  2(l  +  A  )  / 


NOTE:  The  second  (e  )  term  acts  like  a  perturbing  term  which  always  detracts 
from  the  great  circle  route. 

Now,,  .assuming  the  positive  solution: 


1  £ 

~2  \  1/2  "  2 


Fr 1 


(i*i2)(F0 


-1  \  0^ 

g(A)  *  sec  j  -  —  q(A)  (q(A)  a  perturbing  term) 


q(«  ■  [-. - 

J  ( 1  *  x  )  /  r 


U  =  1  +  A 
du  a  2AdA 

A2  ■  p  -  1 

-i} 

e  b 


2Jv 


u(u  -  a) 


q(A)  *  |— T7T  tan 
£  (a)1/2 


(a  =  1  +  82  ) 


q(x)  =  rrrinvr tan_1  ( ftp- 
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^  2 

g(A)  =  9Q  ±  <  sec  (j)-  ^  +  g2  yj 


A2  -  B2 

1  +  e2 


NOTE:  g  = 


A  —  -  1 

U2 


2  ,2  \ 

„  e  A  \ 

1'T'~  ' 


So  for  the  first  three  equations: 


/  2  ,2 
P=b  1  +|--i-T 
l  Z  1  +  X 


(A  >  8). 


4>  =  tan  1  A 

a  a  x  I  "I  A  e2 
9  =  0o  t  Jsec 


ft  -i  /  A2  -  S2  \  1/2I 

- - 6jtanl  X - 1_\  (X>6). 

( i  +  r ) 1/2  \  i  +  8  /  J 


For  the  sphere  the  nodal  crossings: 


eM.  .  s  =  e„  ±  «  . 

For  the  slightly  perturbed  spheroid 

«  -9  +  I  /  1  -  ^  ^  .  | 

(lim  A  +  -)  “  °0  “  2  ^  2  (  1  +  32)  i/2  J  * 

The  new  curve  is  slightly  crunched  up,  presumably  to  avoid  the  eouatorial 
budge. 

The  distance  equation  becomes : 


'"/l1 


2  1  +  A2  j  |(l  +  A2  )2  1  +  a2\a2  (a2  -  e2) 


2  A' 


1  +  A 


.  f _ 1_  .  1  +  ’el  A2 

Jm!  (x2-s2)I/2l  2  i*x2 

•  |x2  -  B2  *  B2  *  82X2  -  e2xV]1/2  d 


-  e2 


-ilLf 

1  +  A2  /_ 
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s  =  b 


'h  A* ' 


dX 


(x2  -  s2 


s  =  b  ( 1  +  ( 1  -  e2  )  $2  ) 


2  N  „2  \  1/2  [  ^  +  e^  X2 


XdX 


2  1  *x2/(i  +  x2)(x2  -  e2  )l/2 


let, 

b"  =  b  ( 1  +  (  1  -  e2 )  32  ) 
XdX 


2\  ,2  \  1/2 


s  =  b 


rI(  1  * 


x2)(  x2  -  e2  )‘/2 

The  first  integral  solves  as 
1 


.  b- 1 '.Jf 


XdX 


xMi*x2)(x2-b2),/2  ' 


,2  V  1/2 


tan 


l/x2  -  32  \  1/2 


(l  +  X2)(x2  -  32)172  * 


( 1  +  e2 )  *'  “  v  1  +  e 

For  the  second  integral 
..  h-  S2  f  X2  XdX 

rj— 

Let 

M  -  1  +  X2 
du  *  2XdX 

_  u'  e2  f  u  -  1  du 


2  r 


du 


i(  W  -  a  ) 


I7T 


-  b'  1 -  f- 

u 


du 

2?  \ l/2 

(  y  -  a  ) 


(1  *82)1/2 


tan 


•i/x2  -  e2\1/2 


-  b 


1  +  3 

2  „2  \  i/2'1 


T" 


b.-.f.l.17,1,  1  tan-L(>-^ 

a“  T72  l  1  .  s2 
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Combining  terms  yields: 


As 


(l  *62)1/2 


,  ,  .2  .2  \  1/2 
/tan-1 1  —4--  ]  I  1  ■  8 


i  +  g 

2  ,(x2-g2)  1/21  x' 


2 

+  T- 


4(l  +  g2  )  1/2  T  +  X2 

b"  =  b  (l  +  (l  -  e2)  g2  )1/2  . 

Now 

'  --2  g2  \1/2 


4(l  +  g2); 
(X  >  g) 


(  1  +  g2  )  1/2  \  1  +  g* 

2 

so  that  (up  to  e  terms) 

.2  „2 


=  b  (  1  -  e 


As  =  b I  1  - 


e  g 


1  +  g2 


1  + 


4(1  +  g2  )  1/2  1  +  X' 

This  can  be  reduced  to: 


2  2  \ 
e  _  e _ \ 

4(1  +  g2)  / 

} 


tan 


"(fr/n 


(x2-62)1/21  . 


As  =  b  <  1  + 


tan 


—  1  /  X 2  -  g2  \l/2 


4 ( 1  +  g2  )/  \  1  +  g 


4(l  +  g2  )l/2 


■  U2-S2)  l/21|Ai 
2 


1  +  X* 


Now  turn  to  the  task  of  finding  S(X)  in  order  to  obtain  instantaneous  heading 
information. 

With, 

pB,‘h(1^TT?) 


<j»(X)  =  tan  1  X 


9(X)  «  0Q  +  sec"1  j  -  |- 


tan 


-l/x2  -  S2\l/2 


(l  +  82  ) 122 
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6(x)  "It L  * p  ar**  * p  sin  *ar  5e 


d<J> 


d6  * 


Eventually  the  ratio  of  two  components  must  be  evaluated.  The  vector  l?(X)  is 

2 

divided  through  by  p  and  expanded  (up  to  e  only) : 
fi(X)  s 


2  X 

1  4, 

1 

1 

Le  ( i  >  x2 ) 2  j 

XP 

J  +  *2_ 

19  + 

(i  TTWJl !-  ij1/2 

1  -  £_ 

2 


1  +  X' 


Since  terms  up  to  e  must  be  retained,  note  that  the  ip  term  does  not  vanish. 
A  look  at  the  geometry  shows  immediately  what  must  be  done  next. 


N0TE;:  i,  does  not  point  along  the  surface,  but  into  the  surface. 

9 

Next,  a  vector  S  (X)  is  needed,  this  vector  is  normally  measured  on  the  sur¬ 
face  (ie,  t)(X)  •  xp  =  Q).  What  is  needed  is  a  rotation  about  the  ig  axis, 
the  angle  is  determined  by  setting  '5  *  ip  =  0,  so: 


5'  (x) 


cos  £  sin  £  0 


1  +  \ 

-sin  £  cos  £  o|  I  e2X  ( 1  +  X2  )  2 


1  March  1988 

First  set: 

0  •=  -sin  £ 


TECHNICAL  NOTE 


(  V) 

+  cos  K 

\  1  +  x2  J 

1  \ 

e2X 


(i  -  X2)2 


-1  2 

tan  £  =  — ; — —  (0  <  X  <  »  ,  ?  <  tan  (e  )  ,5  bounded  and  small). 


1  +  X' 


Now,  since  5  is  very  small,  let  cos  $  ~  1,  sin  £  ~  tan  £ 


e2X 


1  +  X 


fi'(X)  =  - —^r  +  ih  1  eh 


1  ♦  X2  1  .  X2  \(l  *  \TW)^  *  (l  .  A2’)  1/2  '(X2  -  62)I/2 


yf:n  -  1  >  _ _  /,  e2  X2 

'  1  ♦  X2  +  (l  «■  X2  )  1/2  (x2  -  g2  )  1/2  \  2  1  .  X2j  0 

(eliminating  e  terms). 

It  would  appear  from  this  analysis  that  justification  exists  for  keeping  the 
old  t)(X)  and  just  dropping  the  heading  term. 

Now  consider  the  heading  angle  <J».  As  usual,  the  work  is  performed  in  quadrant: 

I. 


tan  =  - -  =  1  - 


2  2 


(X2  +  l)1/2 


0  •  i 


$ 


1+X2/  (x2-32) 


2 „ 2  \ 1/2  * 
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(3)  Slightly  oblate  spheroid  solutions. 


p(X)  «  b  1  + 


e2  X2 


r777 


<|>  (X )  =  tan  1  x 
0  (X )  =0 


As  (X )  =  b 


-1  ■ 


(x2  -  s2)l/2' 


«  Al 


4(1  +  e2)/  4(1  +  e2)  x2  +  1 


4<(X)  =  tan 


I.  e2  A2 

i1  ^TTT 


A2  ♦ 1  \1/2 
A2  -  8  / 


(A  >  8). 


The  most  dramatic  difference  between  spherical  and  oblate  spheroid  geodesics 
can  be  seen  by  comparing  equator  to  equator  trajectories  in  the  oblate  tran¬ 
sition  region. 


Consider  (northern  hemisphere  reference  only) 


For  equator  to  equator  trajectories  where  A0  <  51  ^1  -  jthe  spherical 
geodesic  is  the  same  as  the  oblate  geodesic. 
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For  equator  to  equator  trajectories  where  ir  >  A8  >  ir  J  ,  the  spherical 

geodesic  remains  along  the  equator,  whereas  the  oblate  geodesic  takes  a  route 
similar  to  that  of  the  dotted  line  [actually  there  are  two  geodesic  routes 
(north  and  south)  by  surface  symmetry  about  the  equator]. 

For  equator  to  equator  trajectories  with  A0  =  ir  there  are  an  inifinite  number 
of  spherical  geodesics,  there  are  only  two  oblate  geodesic  routes  (one  over 

each  pole).  In  general,  when  pQ  =  -p  ,  for  the  spherical  case  there  are 

infinite  geodesics,  on  the  oblate  surface  there  may  be  only  two  geodesics 
(over  the  poles),. 

Thus,  another  effect  of  oblateness  is  to  remove  the  singularity  in  the  number 
of  geodesics  which  can  occur  (exception:  for  the  oblate  case  there  is  at 
least  one  pair  of  points  which  produce  infinite  geodesics,  north  pole  to  south 
pole). 


(4)  Comparison  of  FROST  to  integrated  great  circle  route.  Spherical 
minima  routing  gave  a  relation  between  <j>  and  0,  such  that  the  route  was 
confined  to  a  plane.  If  this  route  is  integrated  on  the  slightly  oblate 
spheroid,  the  following  equations  are  obtained: 

P(x)  -  b  (i  +  £. \ 

\  2  1  +  X2  / 

-MX)  »  tan  1  X 

0TX)  =  0Q  ±  sec  1  j  . 


Then, 


and, 

2 

*Up  to  e  only. 
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i|>*  =  tan  1  0 


X2  +  1 

x2  -  e2, 


1/2 


Referring  to  the  geodesy  example  in  Appendix  I  (New  York -Leningrad) ,  the 
integrated  great  circle  parameters  are: 


0  «  0.46113733  . 

X  =  1.130294  . 
a 

Xb  =  0.577350  . 

The  appropriate  level  of  significance  is  ~  e  . 

The  FROST  distance  measurement  gave  AS  s  6735.37  km. 

The  integrated  great  circle  measurement  gives  AS  =  6735.37  km. 

This  result  is  not  surprising  if  the  difference  between  the  two  routes  is 
considered : 


In  both  cases  the  value  of  &  gives  the  minimum  colatitude  the  route  will  take. 
For  FROST 3  24.7379°.  For  integrated  great  circle  =  24.7562°. 

The  difference  between  the  two  comes  out  to  be  approximately  1.1  nautical 
mixes * 


*Up  to  e2 


only. 
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Does  this  result  challenge  the  validity  of  FROST?  No,  rather  it  certifies  it. 
This  particular  example  has  starting  and  ending  points  far  north  of  the 
equatorial  bulge  (the  whole  route  is  far  from  the  equator).  It  is  in  cases 
like  this,  where  the  values  of  X  and  0  are  relatively  small  that  one  would 
expect  differences  between  the  two  routes  to  be  vanishingly  small.  On  the 
other  hand,  routes  predominantly  close  to  the  equator  should  show  greater 
distance  and  routing  deviation  (consider  two  points  on  the  equator  179.9° 
apart  in  longitude) . 


K.  Carlton-Wippern,  PhD. 
Mathematician 
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10.0  PROJECT  FROST 

10.1  The  need  for  geodesic  and  surface  navigation  algorithms 

There  exists  a  continuing  military  requirement  to  review  and  upgrade  (as 
needed)  the  algorithms  used  to  perform  geodesic  measurement  and  perform  sur¬ 
face  navigation.  These  algorithms  affect  the  manner  in  which: 

a.  Mapping  agencies  measure  distances  between  points  (along  a  geodesic 
or  curve  of  minimum  distance) . 

b.  Long  range  aircraft  (LRA)  navigate  on  the  reference  surface. 

c.  Long  range  naval  operations  (surface/submarine/ air)  perform  surface 
navigation. 

d.  Long  range  cruise  missile  surface  guidance  is  performed  during  iner¬ 
tial  flight,  phase. 

e.  Strategic  defense  initiative  surface  interceptor  problems  can  be 
addressed. 

(1)  Airborne  platform  positioning  and  navigation. 

(2)  Reentry  intercept. 

(3)  Interception  of  compressed  trajectory  vehicles  via  kinetic  kill 

vehicles. 

f.  Other  (NASA  type  applications). 

10.2  Current  methodology 

Use  spherical  trigonometry,  with  a  compensated  equal-area  radius, 
a.  Basic  spherical  trigonometry  relations. 
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(1)  Law  of  tangents. 


cos  Avfi  =«  cos  $Q  cos  ^  +  sin  <J>0  sin  <*> ^ 
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(3)  Law  of  sines i  Law  of  sines  generally  used  for  heading/ 
directional  information. 


C 

b.  Current  method.  Uses  spherical  relationships  and  equal-area  radius 
value. 


(1)  The  equivalent  area  sphere  to  oblate  spheriod  relation  is  as 

follows: 


4wl2 


2ira 


An 


1  +  e 

1  -  e 


Rearranging : 

a2  -  §  ♦  Je  (  In  M  +  s|  -  In  |1  -  e|  ) 

2  4*[a2+  7  (-e  -  4  +  -I  +  e  +  £  (assuming  e2  small) 
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(2)  The  above  approximation  for  the  radius  value  is  to  be  used  with 

2 

Naperian'  trigonometry  (assumed  to  be  valid  up  to  e  correction). 

10.3  Problems  with  the  current  methodology 

a.  On  an  oblate  spheriod  the  radius  vector  magnitude  varies  from 
3  <  I  r  I  <  a. 


b.  Now  since  e,  the  eccentricity  can  be  defined  as  e  *  (  1  - 
then  rearranging  ' 


$ 


a  ■  b(1  -  e2fl/2 


and  if  e  «  1  then, 
2 


a  J  b 


(-!)• 


therefore: 
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—  »r 

b<  |?|<  b  (i  ♦  §  )  . 

c.  For  the  equal  area  compensation 


This  value  of  |  r  |  is  assumed  to  be  used  over  the  whole  spheriod. 

d.  A  problem  arises  in  the  level  of  accuracy  required  for  precision 
measurements.  The  so  called  equal  area  method  should  have  an  accuracy  limita- 

tion  on  the  order  of  (or  no  less  than  an  order  of  magnitude  less  than)  e 
2 

(on  the  earth,  e  *  .0067).  For  quarter  turn  measurements  (the  distance  swept 
out  by  a  90°  central  angle,  this  to  be  redefined  later),  errors  from  -5  to  +11 
kilometers  can  be  obtained. 

10.4  Introduction  to  FROST 


a.  Start  by  rederiving  the  spherical  case  geodesics,  considering  each 
necessary  component  separately,  and  relating  each  component  via  a  parametric 
variable,  A. 

(1)  For  the  spherical  case: 


P  (A  )  *  a 


9  (A )  =  tan  1  A 


0(A)  «  0Q  ±  sec  1  j 


9(A)  »  tan 


-1 


A2  +  1 

A2  *  32 


1/2 


(synthetic  law  of  sines) 


As  (A)  -  a  tan 


■l/x2  -  e2\1/2 11 

l  1  *  82  / 


(synthetic  law  of  cosines). 


(a)  It  can  be  shown  that  upon  elimination  of  A,  the  above  rela¬ 
tions  reduce  to  standard  spherical  trigonometry. 
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(b)  Spherical  case  quadrant  analysis  for  1J1 
QUADRANT  I  QUADRANT  II 

i 


DIRECT 


RETROGRADE 


for  $  the  heading  or  direction  angle,  tan-  1  0 
the  following  terms: 


X2  +  1 
X2  -  3 2 


1/2 


=  q ( X )  produces 


Quadrant  I 

..Quadrant  II 

q(X) 

ir  -  q(A)  DIRECT 

it  .  +  Q  (X  ) 

2tr*  -  q(X)  RETROGRADE 

(2)  For  the  slightly  oblate  spherical  case  assume  the  eccentricity 
value  e  is  small.  Then  using  calculus  of  variations  and  solving  for  the 

2  4 

geodesic  curves,  keeping  terms  up  to  e  only  (eliminating  e  and  higher  terms) 
obtain  the  following  relationships: 

p  (X )  *  b  [  1  +  -I  •  — - — ~  \ 


1  +  X‘ 


$  (X )  =*  tan  1  x 


9(X)  -  0Q  ±  sec"1  i  -  # 


6  2  (l  ♦  02)1/2 


tan 


-1  /X2  -  82  \1/2 

,1  ♦  02 
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♦  (X)  =  tan 


2 

_e 

2 


_i!_\ 

1  +  X2/ 


8 


As  (X) 


2 

o 

4(1  ♦  e2 ) 


tan 


4(1  .82)  1/2 


(x2  -  82  )  1/2 


X  +  1 


(a)  <j>  in  spherical  coordinates  is  defined  as  the  geocentric 
colatitude  variable.  The  geodetic  latitude  variable  is  related  to  X  by 

L(X )  =  cot-1  [(1  -  e2)  x]  . 

(b)  In  the  above  equations,  in  the  limit  as  e  +  o,  the  spheri¬ 
cal  case  solutions  are  recuperated. 

b.  FROST  makes  a  first  order  attempt  at  solving  for  the  geodesics  on  a 
slightly  oblate  surface. 


(1)  It  assumes  the  surface  eccentricity  (e)  is  so  small  that  terms 

4 

of  the  order  of  e  and  higher  may  be  discarded.  This  seems  a  reasonable 
assumption  since: 


turn. 


(a) 


2 


e 


corrections  are  from  -5  to  +11  kilometers  per  quarter 


4 

(b)  e  corrections  are  on  the  order  of  tens  of  meters  per 
quarter  turn,  assumed  to  be  well  within  the  noise  of  the  actual  surface  fluc¬ 
tuations  of  the  geoid;  therefore, 

4 

(c)  Higher  order  terms  (greater  than  e  )  can  be  absolutely 
discarded  for  the  earth  for  any  reasonable  accuracy. 


(2)  The  above  assumptions  form  the  basis  for  the  limitations  on  the 
use  of  FROST. 


10.5  Comparison  study 

a.  In  comparison  with  the  spherical  relations  note  an  apparent  crunching 
up  of  the  curve  for  9  at  the  nodal  crossings. 

9(lim  X  00 )  *  9q  ±  j  sphere 

d 

1  +  62) 


■jyjj oblate  correction 


0  (lim  X 


9  ±  — 

o  2 


2 

e 

T~ 


also,  iji(node)  yields: 
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♦  dim  A  *  •)  =  tan  8  sphere 
-1 


♦  dim  X  »)  =  tan 


i'-ty 


oblate. 


NOtE:  ♦  dim,#  ®)  =  in  both  cases. 


As  |  =  a  •  rr  sphere 


A=8 


ls  1  x=e  ■  b 


(i*(i  -_£wyii,  i 

(  1  *  SZ  ) 1/2  X  2  4(1  * 


82) 


ir_ 

2 


A  1  x=®  „  irb  /  4 

As  ,  _a  S  t-  1  ~  o 


2  „2 

e  0 


X=8  "  2 


1  + 


e 


2  1  ♦  82A  7  4(1  +  62  )., 
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oblate. 


s  ^  (l  +  f-  •  — - — -  )  (keeping  up  to  e2  terms  only). 

2  \  1  +  e2  / 

(1)  At  first  glance,  this  may  look  like  nonsense,  as  8  gets  larger, 
the  distance  traveled  is  smaller,  but  look  at  9(A).  As  8  gets  larger  A9  =  ir 
is  not  transverse  so  that  this  result  does  make  sense. 

(2)  Parametric  geodesics,  synmetric  about  an  arbitrary  9q. 

(a)  Sphere. 
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(b)  Slightly  oblate  spheroid: 


b.  Comparison  of  distance  formulas.  Oblate  spheroid,  perturbed  spheri¬ 
cal  routing  with  great  circle,  equal  area  routing.  Comparison  will  be  from 
the  high  point  in  the  trajectory  (or  maximum  latitude)  (  X  ~  3)  to  the  equator 
(X  *  •)  so  that  A  is  evaluated  as  8  <  X  <  •. 

5 


(1)  ^or  the  perturbed  solution. 


A0  = 


1  "  ri 


4sn  -?(1 


,(l  *  S2  )  1/2 

o2  \  1/2 


1  8' 


1  + 


s  . 

7" 


4(l  +  02  )i 


X*8 
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(2)  For  the  spheroid  approximation: 


As  ,  TT  ...  IT  AQ 

cos  —  =  cos  <j>  cos  ■=■  +  sin  <j>  sin  ■=■  cos  A0 

3  4  4 


As  =  a  cos  1  (sin  <J>  cos  A8] 


but  cos  A0  *  sin  . 


(777 J57?  • 


It  is  convenient  to  call  B  in  these  cases  tan  i  ,  the  minimum  colatitune 

o  0 

value  for  a  given  geodesic,  so  that: 

Asa  *  a  cos  1  sin  <j>0  •  sin  ^  sin  <j>Qj  s  lr  ^1  -  f-  sin2  <J>0  j  . 

*  .  ir  (  ,  2  .  2  \  1/2  ( ,  e2  e2  2  ,  \  „  nb  ( .  .  e2  2  ,  \ 

Asn  3  b  1  (  1  "  e  sin  ^oj  [1  T“T"c°s  ^o)=7"  |1  T"  cos  hj 


2 

(up  to  e  terms). 

Now,  with  e2  =  .0067227, 

a  *  6371.227  km 
b  *  6356.912  km  . 
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1968 

*r 

(3)  For  the  quarter  turn  (defined  as  8  < 

Asn 

Asa 

0 

10002.196 

10007.899 

10 

10001.690 

10006.885 

20 

10000.233 

10003.964 

30 

9998.000 

9999.489 

45 

9993.805 

9991.079 

60 

9989.609 

9982.669 

75 

9986.538 

9976.513 

90 

9985.414 

9974.259 

All  figures  in 

kilometers. 

As  =  As 
n 

a  at  *0 

=  35.568° 

Max  difference 

As  -  As 
n  a 

positive  »  11.155  km 

As  -  As 
n  a 

negative  =  5.703  km 

The  approximation  is  good  to  two  significant  digits 

c.  The  comparison  of  FROST  versus  spherical  trigonometry  is  shown  in 
figure  1. 

10.6  Geodesic  measurement  example 

a.  Measure  distance  between  two  theoretical  bases: 

(1)  Base  A:  4125000  N,  7321146  W  (New  York)  geocentric  latitude. 

(2)  Base  B:  6020000  N,  2928000  E  (Leningrad)  geocentric  latitude. 

(3)  Alongitude  -  10229146 

(a)  4126918  N  geodetic  latitude. 

(b)  6021670  N  geodetic  latitude. 
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Figure  10.1  Measurement  comparisun  for  a  quarter  turn 
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b.  Comparison  of  distances  measured. 

2 

{  CD  Spherical  trigonometry,  equal  area  e  compensation,  geocentric. 

-  .As,  =  6744.0840  km. 

2 

(2)  Spherical  trigonometry,  equal  area  e  compensation,  geodetic. 

As  =  6716.2297  km. 

(3) .  FROST  (6  =  .460749734,  ,0O  =  722318). 

As  =  6735.3776,.  km. 

(a)  Spherical  trigonometry  w/e  compensation,  geocentric. 

Let> 

A  be  at  longitude  o° 

B  at  longitude  10229146 
A  *  {.748955721,  0,  .662620048} 

B  =  {-.111749248,  .487352137,  .866025404} 

0  =  60264952234 

As  =  b  (l  +  )  (1.058533854) 

=  6744.0840  km. 

NOTE:  0  =  cos'1  (A  •  B) 

2 

(b)  Spherical  trigonometry  w/e  compensation,  geodetic. 

A  *  {.746732702,  0,  .665124253} 

8  *  {-.111184443,  .484888955,  .867479522} 

0  =  60239902851 

As  ■  b  (l  +  |-)  (1.054161913) 

=  6716.2297  km. 
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NOTE:  9  -  cos-1  [t  •  §) 


(c)  FROST  calculation  (9Q  *  7J2318,.  0  *  .460749734) 


A,  Xft  =  1.130294,  XjJ  =  1.277565 
8,  XB  =  0.577350,  Xg  =  0.333333 

.  ,  -1  /1. 277565  -  0.2122903  l,1f 

Asa  *  ta  tan  { - 172122903 - )  V 


0.006722 


0.006722 


4  (1.2122903  )1/Z 


L 1.277565  -  .2122903 
2.277565 


.  ,  -1  f. 333333  -  0.2122903  \  l,1f „  .  0.0 

Asb 3  b  tan  v - 1 .2122903 — )  V  +  4 rr. 


006722 

.2122903] 


'j  1/2 


0.006722  (  .33333  -  .2122903  ) 

'  4  (1.2122903)^'  - ™ 

As  -  b  {.754145026  -  .000691664}  -  b(. 753453361 ) 

Asb  -  b  {.306480974  -  .0003982597}  *  b( .306082714) 

As j  -  b(1. 059536075),  with  b  -  6356.912 

«  6735.3776  . 

Given, 

p  *  .460749734 

X  -  1.130294 
3 


Xb  -  0. 577350 


calculate  A9^.: 


A9  ■  sec  *  g-  - 


•'/»s  B2V/2 


(1  -S2)1'2  \  1  +  8 
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A0L  =  65.94351463  -  .06060861  2  65.002826 
A9r  *  37.05647301  -  .024663507  *  37.031809 
A8t  =  102.9146(35)  . 
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APPENDIX  II 

20.0-  THE  TORUS 

The  torus  is  a  very  complicated,  multiple  connected  surface,  and,  as  a  result, 
one  should  not  expect  nice  closed  form  solutions. 

What  is  presented  here  is  just  an  outline  of  how  one  might  go  about  solving 
for  the  geodesics  on  such  a  surface.  The  solutions  presented  are  strictly 
formal. 

Oefine  the  surface  model: 


(r  -  a)2  +  z2  =  b2  (a  >  b) . 
Then  let. 


This  can  be  simplified  somewhat  by  letting  \  *  tan  then 
r  *  a  ±  b  cos 
2  =■  ±  b  sin  <J> 
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8  (#)  =0  +  6 e  f _ ^ _ 

°  J  r(  r2  -  c2  )A  1 

4SW,‘b/(r2-!2)^d* 

S(<|0  =  ±  iz  ±  tan  ± 


1  ± 


ce 


(f2- 


e  cos  ijf  ) 
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APPENDIX  III 

30.0  A  STUDY  ON  SURFACE  CURVATURE 

Surface  curvature  may  be  evaluated  using  the  parametric  geodesics  and  the 
appropriate  Frenet-Serrat  equations  because  the  geodesics  for  a  given  surface 
contain  within  them  qualitative  as  well  as  quantitative  information  regarding 
the  surface  itself.  The  parametric  geodesics  form  a  natural  network  or  grid 
which  define  the  surface. 

This  evaluation  is  started  by  presenting  the  Frenet-Serrat  equation  for  cur¬ 
vature,  evaluating  the  Christoff el  symbol  of  the  second  kind,  showing  the 
parametric  evaluation  for  the  symbol,  then  proceeding  directly  to  the  spheri¬ 
cal  case. 

i  .  d  i 

If  x1  are  the  coordinate  variables,  X1  »  -p  the  unit  tangent,  then 


dX*  r1  a  dx®  i 

XT*  +  *  *  Ky-1  (where  K  is  the  curvature,  p1  is  a  unit  vector), 

ds  |a>f}  ds 


-r  , 


*  -r  sin"  <{>  , 


C-»- 


r2 

*  -sin  $  cos  $  , 

*3 ,3 


cot  $  .  All  other  f1  =0  . 


For  the  sphere  (r  *  a  )  and  , 


-a  , 
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r 

•2,3 


1 

X 


30.1  '  Spherical  case 
Letting 


2  -l 

x4  -  tan  X 


x3  ■  9Q  +  sec  1  j 


(choosing  the  positive  side) . 


.  O  .1  \ 

NOTE:  The  result  is  the  same  if  x  «  8n  -  sec  £  is  choosen. 

o  p 

Then, 

X  ^  •  0 

,2  1  dX 


1  +  X 
3  1 


2  3s 


.  dX 


X  —  -  1 

U2 


1/2  ds 


and  with 
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l 

a  * 


Evaluating, 


dX2  d 

F  3T 


IO 

1 

ra 

fo 

)1'*  \ 

a 

Cl  +  02  ) 

l/rW 

dX 

3? 


1  March  1988 


TECHNICAL  NOTE 


=  0  . 

dX3 

Evaluating  is  somewhat  more  difficult.  Remember  that  the  great  circle 
route  in  parametric  curvilinear  coordinates  is  a  multiple  defined  curve 
9  *  9q  ±  sec  1  j  ,  based  on  a  multivalved  function. 

Consider: 


The  first  derivative  is  ±  ^  ,  depending  on  which  side  of  the  curve  (or, 
batter  said,  which  branch)  is  chosen  for  usej  however,  the  second  derivative 
d29  , 


d\ 


is  both  equal  in  magnitude  and  direction  (all  concave  down)  on  either 

d20 

side,  so  that  a  straight  forward  brute  force  evaluation  gives  contribu- 

dX2 

tions  from  both  branches.  The  term  must  be  divided  by  two  for  use  in  the  cur¬ 
vature  equation. 


Continuing  with. 
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is  known  to  be  generally  true,  not  only  for  the  great  circle  routes  but  for  the 
surface  of  a  sphere  [constant  curvature,  independent  of  where  the  observer  is 
(0O,  X)  or  in  what  direction  the  observer  is  looking  (0)].  This  is  an  impor¬ 


tant  result.  It  also  points  out  what  would  have  happened  had  not  1/2  of  the 

dA3  3 

term  been  taken  when  evaluating  the  curvature  (Ky  ?£  0  but  a  function  of 

X,  3  which  would  not  make  sense). 

For  completeness  lets  turn  our  attention  to  evaluating  the  torsion  of  great 
circle  routes. 

Fortunately  it  is  a  trivial  matter  to  show  that  t,  the  torsion  for  the  para¬ 
metric  spherical  geodesic  is  x  =0. 


From 


rv 


now 
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ds 

and 


What  is  needed  to  be  shown  is: 

f  yaX6  +  KX1  -  0  =  xy1  . 

'ag 

NOTE:  a  =  1  always,  otherwise  y 
For  i  ■  1, 


=  0. 


r. 


Tyl  »  0  . 
For  i  *  2, 

r2-i 

a 

1,2 

so  that 


(aX'1)*2  *  aX*  '  0 

Ty2  ■  0  , 

For  i  -  3, 


1.3 


TECHNICAL  NOTE 


112 


1  March  1988 


TECHNICAL  NOTE 


tu3  =  Q 
QEO. 

For  the  minima  solutions  on  the  sphere  there  is  zero  torsion  and  constant  cur¬ 
vature,  independent  of  $  or  X.  It  is  probably  primarily  for  this  reason  that 
the  sphere  is  so  clean  to  work  on. 

The  end  of  this  section  will  concern  itself  with  the  evaluation  of  curvature 
on  a  slightly  oblate  planet.  Coordinates  to  be  used  are,  of  course,  spheri¬ 
cal. 
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K|i 


.  (3g2  +  2xV  -  2X2  ) 


b2  (  t ■  ♦  X2  )  X  ( 1  +  B2  ) 


■r  =  o  . 

Combining  all  terms,  the  curvature  can  be  written  as: 


K  u1  ♦  K,  u1 
(TO  Vl 


where  is  an  e2  perturbing  term. 


'ri 


lRl  -(sij(Kui)(K“J))1/2)  • 


30.2  Conic  curvature 
P  *  aX 


!j>  =  sin  1  a 


a  a  +  1  „  "1  *<* 
0  -  90  1  -  sec  — 


^the  positive  branch  is  chosen  and  as  in 

j2g  \ 

the  sphere  take  one-half  the  value  of  j  . 


dp_ 

dX 


£  ■« 


d6  1  1  c 

r^r—vm 
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the  case  of 
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Ku=0  +  0-  a(l  -a) 


2  n1/2 


2  4.4  * 

a  a  X 


Km 


S  .  .  c(«2X2  -c2)1/2  ,  1  /(„2X2  -  c2  ) 1/2 

2  3.4  +  aT  1  aX 

a  a  X 


„/  „2,2  2  N  1/2  /  : 

.  -  cta  x  ~c  )  +  pA«. 


2,2 

aa  X  , 


+  0 


2  3,4 
a  a  X 


2,2  2  \  1/2 

X  -  c  ) 

2  3,4 
a  a  X 


=  0. 

So  that  for  the  cone. 


Km1  =  c2«  ( 1  -  «2  )  1/2 
^  2  4  4 

a  a  X 


0 

-1  |  (lim  a  1,  K  =  0) . 

0 


Evaluating 


|Km 


i  |  ,  1 .  c2a  ( 1  -  a2  ) 

•  a  U  1 


1/2 


•V 


NOTE:  |K|  -  ~  . 


It  is  easy  to  show  that  the  curvature  f  a  helix  (the  parametric  cylindrical 
geodesic-)  is  in  the  Ir  direction  only.. 

For  the  three  cases  (sphere,  cone  and  cylinder)  where  it  was  shown  that  the 
surface  constraint  was  written  of  the  form  u^  *  constant  (i  «  1,  2  or  3>  u,  a 

curvilinear  variable),  the  curvature  Km1  of  the  represented  sur-ace  ' as:;  rin* 
it  existed)  had  a  component  ONLY  in  the  &.  direction. 


Whether  or  not  this  is  a  general  property  is  beyond  the  scope  of  this  tech¬ 
nical  note. 
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APPENDIX  IV 

40.0  SURFACES  OF  REVOLUTION  AND  LONGITUDINAL  ASYMMETRY 

Calculation  of  geodesics  on  surfaces  of  revolution  and  first-order  pertur¬ 
bation  due  to  longitudinal  asyrrmetry  (analysis  in  spherical  coordinates)  can 
be  accomplished  as  follows. 

Starting  with  the  spherical  surface  metric: 

assume  the  surface  can  be  written: 

P  *  p  C4» )  ,  [■  p  (X ) ]  ,  a  surface  of  revolution. 

Then  with-  4>  =  tan  1  X  ,  9  =  g(A)  obtain: 

[('•a«‘KTT7)'-77?  .-f’* 

Let(r If-)* q(x)  th8ns 

As  =  fp  (A )  |~(1  +  q(A) )  •  f - -  \  +  — - — -  g'2  dX 

J  L  V  1  +•  X2  /  1  +  A2 


Since  the  partial  of  the  integrand  with  respect  to  g  »  0,  then,  as  done 
before,  let 

T  /  ,  SZ  ,2  _T-l/2  ,2 

x  I  /  .  \  /  I  \  \  -*2  I  A  ^ 


c  -  Pft)  •  (1  *  qttnf - g'2 

'  1  +  X  /  1+X 

where  c  is  the  first  constant  of  integration. 


1  +  X' 


Solving  for  g  gives: 


1  +  q( A) 
_2  # 


X2 


L  p  u) 


(i  +  x2) 


8  (X )  «8, 


±  c  f _ 3 _ II _ 

Jx‘pu)  Lx2-4(1+x2). 
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0O ,  the  second  constant  of  integration. 
For  the  sphere  p  =  constant, 
q  =  0, 


letting  —  »  5 


8'  -*£•  [x2(l-c2)-c2]1/2  . 


Letting  — - — -  =  02 
1  -  c 


’  Is2  J 


9  ■  9g  ±  sec  1  j  as  before. 


For  a  slightly  oblate  planet,  where  e4  terms  are  significant,  but  not  e6 
+  terms,  proceed  as  follows: 


<KX)  *  tan  1  X 
calculate  8 (X )  by  noting 


q(X)  »  e  • 


1-  -  e' 


1  ♦  X 
,2  .  X 


2  \  2 


1  +  V 


*(s-r 


g  -II  -  a' 


2  X 


2  \  1/2 


1  +  X 


1  +  e 


1  ♦  X' 


1  +  X" /  J 
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H?-1 


1/2 


g  S  1  -  f-  • 


2  ,2 

e  X 


1  ♦  X 


+  e 


2(l  +  X2)2  8  (  1  +  X2  )  2 


9  ■«„  i'jseo"1  82)l72 


XdX 


e  tan->  (i^£  Y/2 


1  +  B 


X3dX 


)1/2-(  1+*2)2  J4(^-i)1/2  •  (1  ♦  x2)2 


This  solves  explicitly  as : 


*»>■'.  r  -  3- 


tan 


-l/x2  -  82  V/2 
1  +  a2  / 


58 


X2  -  a2  u/2 


8  2 ( 1  +  a2) 1/2  \  1  +  x2  \  1  +  a2  j  (1  +  a2) 


tan 


■i/a2  -  S2\1/2 


tan 


-i/a2-b2\,/2 


1  +  82  /  /  (l-B2)1'2  \  1  +  8 


ip(X)  and  As (X )  can  now  be  solved  for  by  the  usual  method: 


-1  J  /  x2  +  1  \1/2 

*(X)  -  tan1  a  •' 


is  »  b  1  •  s 


2  8 


x2  -  a2 


1/2 


1  +  8 


tan 


2  .2 

,  e  X 

1  ~  2  ’  ,2 

1  +  X 


•i/a2  -  62\1/z 


1  *  8 


4(l  »B2)l/2  ’  1  *  A2 


■  (*2  -  B2  ) 1/2  ♦  2sl  .  tan-2  (  ^  -  92  Y/2 

8  1  *  82  I 


e4 

1 

8 

_1  +  X2 

A2  -8 2\1/2  .  1 


tan 


-l  /  X2  -  a2  V/2~ 


( 1  +  a2 )  \  1  +  a 


e 

64 


A2  -  B2W2 


(l  *  A2  )2  \  1  *  B2  /  (l  ♦  X2)(l  ♦  B2 )  \  1  *  82 


A2  -  B2\I/2 


tan 


-1/x2  -  a2  V/2 


(l  ♦  82  )2  \  1  ♦  8 


J  J  X, 
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Longitudinal  asymmetry  is.  a  very  complex,  highly  nonlinear  problem,  solvable  pri¬ 
marily  via  successive  approximation.* 

Start  by  assuming 

P  *  a  ( 1  +  eg)  where  e  «1 


In  the  analysis  assume  the  primary  interest  is  in  working  on  a  small  surface 
portion,  or  sectoral  slice,  where  the  above  assumption  is  approximately  true; 
also,  assume  working  in  the  low  to  mid  latitude  region  is  of  primary  interest. 


thus  excluding  analysis  in  the  polar  region,  ie. 


e  ~  0  and  e  << 


X2 


for 


all  X  in  the  domain  of  interest. 

On  the  earth  such  an  example  would  be  the  geodesic  from  St  Louis,  Missouri 
(elevation  ~  1000  feet)  to  Denver,  Colorado  (elevation  ~  6000  feet). 

Continuing: 


As  *  a 


The  integrand  can  be  written: 


(1  +  eg) 


1  +  X' 


♦  e2 


'1 


1  +  X' 


1/2 


dX 


which  can  be  further  approximated  as: 

i  1/2 


(1  +  eg) 


— - —  )2  *  _ g'2 

1  +  X2  j  1  +  X2 


dX  . 


•Reference  Einstein,  Law  of  Gravitation,  pgs.  328-334 >  Lass,  Vector  and  Tensor 
Analysis. 
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In  this  example,  letting  the  integrand  be  6,  we  see  immediately  ^  0  , 

Therefore  we  must  use: 

d  36  _  36  _ 

d*  3g'  3g  " 

and  collect  terms  in  e,  setting  up  for  a  successive  integral  approximation. 
Continuing, 


Rewriting  the  equation  as  — t-  - 

8g 

(Since  as  e+0  we  must  obtain  the  unperturbed  spherical  solution,  and  letting 


The  system  becomes; 


Start  by  setting  e  =  0,  solving  for  g  obtain: 


a  .  -1  X 

g0  ■  90  *  380  J  • 


Substituting  this  on  the  right  side  of  the  equation,  obtain: 
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I  .  -IX 

•  <  c  0Q  ±  c  sec  j  + 


X  _/  .  .1/2  . 

- _  g  -  c  [  1  =  e(X) 

1  +  \ 


1  ♦  x2/  (x2  ♦  1  )  /  xf  _  1 


where  e (X )  = 


{oeo±osec'1T*tan'1(7f/-)1/  j  • 

The  next  iteration  for  g  (=>  g^  'can  now  be  set  up  and  performed  numerically 
(note  e(X)  has  a  ±  term,  this  accounts  for  asyrmetry  in  the  iQ  direction). 


Now  obtain: 


,X2(x2  (l  -  c2)  -  c2)g'2  =  c2  +  2c  e (X)  •  (  ]l/2  (e2  ~  0). 

continuing  the  process  of  successive  approximation: 


*2(*2  Cl  -  c2)  -  c2)/2  .  c2  *  2c  .»)  •  -A-  .( 

1  *  X2  \  X  -  3  / 


2  \  1/2 


and  finally: 


g  -  9 


+  C  11  *  e  1  r(X)) 


X  —  -  1 

U2 


where  r(\ )  ■  —  7 — -  —  ? —  (  ~~  —  — 

°(i  *X2  2  U2  -  s2 


c  eo  ten- (Ai^|!)  1/2 


NOTE:  0 


2  c 

~ 


Caution  must  be  used  when  determining  which  combination  of  ±’s  need  to  be 
employed,  and  on  which  side  of  8q.  Once  this  is  accomplished,  <i>(X)  and  As(X) 

can  be  calculated  in  a  straight  forward  manner. 
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PREFACE 

This  Is  the  first  addendum  to  AFSPACECOM  Technical  Note,  Surface  Navigation 
and  Geodesy,  A  Parametric  Approach,  1  March  1988.  This  addendum  will  address 
two  issues  not  covered  in  the  aforementioned  publication. 

First,  a  practical  problem  involving  conic  geodesy  will  be  presented  and  the 
solution  discussed. 

Second,  the  solution  for  the  parametric  geodesic  of  a  paraboloid  (parabola  of 
revolution)  will  be  presented,  but  the  analysis  will  involve  the  use  of  para¬ 
bolic  confocal  coordinates. 

It  is  hoped  this  addendum  will  further  demonstrate  the  value  of  the  use  of 
parametric  analysis  in  solving  problems  of  surface  geodesy. 


THIS  TECHNICAL  NOTE  IS  THE  WORK  OF  OR.  K.  CARLTON-WIPPERN  ANO  DOES  NOT 
NECESSARILY  REPRESENT  THE  OFFICIAL  VIEW  OF  THE  US  AIR  FORCE. 
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I.  AN  APPLICATION  OF  CONIC  GEODESY 

In  many  circumstances  the  surfaces  of  certain  physical  features  can  be 
approximately  modeled  (to  a  very  high  degree)  by  that  of  a  cone.  A  case  in 
point  would  be  Mount  Fuji  in  Japan,  which  will  be  assumed  to  be  in  the  shape 
of  a  perfect  cone.  Consider  the  following  problem: 

A  group  of  Shinto  monks  are  making  a  pilgramage  from  their  home  in  Kofu  to 
their  seaside  temple  at  Numazu.  On  the  route  their  master  has  selected,  they 
will  have  to  pass  across  the  base  of  Mount  Fuji,  from  the  northern  most  point 
to  the  southern  most  point.  The  leader  of  the  noble  group  of  pilgrim  monks 
has  been  charged  with  taking  the  route  about  or  over  Mount  Fuji  which  will 
traverse  the  minimum  distance.  It  can  immediately  be  seen  by  symmetry  that 
there  will  be  at  least  two  routes,  one  about  the  east  face  and  the  other  about 
the  west  face.  Not  withstanding  this  temporary  dilemma  (which  will  be  decided 
by  the  toss  of  a  coin)  the  following  elements  will  have  to  be  calculated/ 
determined: 

a.  The  minima  curve. 

b.  The  total  distance  traveled. 

c.  How  high  up  the  side  of  the  cone  will  the  monks  need  to  travel. 

d.  The  initial  angle  of  attack. 

Given,  Mount  Fuji  has  conic  parameters. 


So  that, 

P  =  (h*  +  r*)1/* 

0 

r 

«  *K 

r  S  r 

(h*  +  r*)*/* 


r 


i 
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Solution.  Reference  paragraph  IV,  Surface  Navigation  and  Geodesy,  A  Parametric 
Approach,  and  obtain  the  parametric  conic  geodesic: 

p  =  aX 

$  =  tan-1  a  (constant) 


-  $  f 

A  diagram  of  the  problem. 
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so  that 


tp.  / Tir  \ 

c  =  “r  cos(r)  * 

this  fixes  the  parametric  curve  (0Q  =  ^  for  west  face  ascent,  for  east  face 
ascent.) 

Now,  p_.  occurs  at  I  8  -  0  I  =  0,  so  that 
mm  1  o 

rXmin 

,  _  min 


From  th’is  obtain 
pm1n3P„cos(?)  ' 

The  distance,  measured  from  ground  level,  traveled  up  the  side  of  the  mountain 
is 

\  -  pmin  =  \  t1  '  005  (?)]  ' 

The  total  distance  to  be  traveled  from  one  side  of  the  mountain  to  the  other 
can  be  evaluated  as  follows: 


1 


So  that,  for  both  ascent  and  descent  parts  of  the  geodesic: 

x/2  x=(P0/a)=Xmax 

X=(pmin/a)=Xmin 


is*2a(x!' 


Now,  since 


i  =  — 
min  T  » 
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Po 

Nnax  “  aT  ’ 


The  initial  angle  of  attack,  defined  as  the  angle  as  measured  off  from  the  top 
of  the  mountain,  measured  in  the  plane  whose  normal  is  parallel  to  the  conic 
surface  normal. 


♦  *  |  0  -  t>  . 

Summarizing  the  solutions: 

Total  distance  traveled  s  2p 

Distance  of  ascent  up  the  mountain  side 

Initial  angle  of  attack  *  |  (1  -  T). 


=  P([i  -  «$(£)] 
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II.  SOLUTION  OF  THE  PARAMETRIC  GEODESIC  FOR  A  PARABOLOID  (SURFACE  OF  A 
PARABOLA  OF  REVOLUTION) 

This  section  will  be  brief,  outlining  only  the  most  necessary  steps  required 
to  get  to  our  objective,  that  being,  defining  a  surface  constraint  and  develop 
the  five  equation  parametric  solution  for  the  geodesic  curves  (the  necessary 
surface  navigation  relations). 

For  simplicity  sake,  the  surface  defined  in  cylindrical  coordinates  is: 

r«-2a(z+§)  , 

where  a  is  an  arbitrary  parameter. 

Parabolic  confocal  coordinates  can  be  described  as  the  set  (£,  n,  0)  such  that 
0  is  the  same  longitudinal  variable  as  described  in  cylindrical  or  spherical 
coordinates, *then  using  p  from  spherical  coordinates  and  Cartesian  z, 

l  =  p  -  z 


T)  =  p  +  Z 

and  the  surface  constraint  equation  can  be  transformed  to: 
r2  =  a(2z  +  a) 
z2  +  r2  =  z2  +  2az  +  a2, 
z2  +  r2  =  p2 

P2  =  (z  +  a)2  * 


p  =  z  +  a 
p  -  z  =  a  =  $ 

hence,  the  surface  constraint  is  equivalent  to  setting  one  of  the  curvilinear 
variables  to  a  constant,  namely  £  s  a. 


Since  the  Jacobian  matrix  is  formed  to  show  curvil inearity  and  obtain  the 
scaling  factors,  it  is  convenient  to  work  with  spherical  coordinates  as  an 
intermedium: 


From 


(p  sin  i  cos  0\ 

p  sin  f  sin  0 

p  cos  *  I 


/  p  sin  i  cos 
p  sin  ^  sin 
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p  -  z  =  5 
p  +  z  =  J| 

2z  =  T)  -  $ 
so  that 

z  =  i(n  -  K)  . 
Now 

p  =  Kn  +  5) 
so  that 


cos  <*  =  -  = 

p  n  +  5 


*  (n*  +  2n{  +  {'  -  n*  ♦  2n£  ♦  5!)‘/> 

-  2(n£)1/1 
n  +  $ 

x  *  p  sin  <t>  cos  8 

*  *<"  +  5>2  cos  6 

=  (tiO1/*  COS  0 

A  similar  relation  exists  for  y,  so  that 
/  x  \  /  (nO1/2  cos  8  \ 

(  y  )  “  (  (^S)1/2  s1’n  9  j 
\z/  \Hn-5)  / 

The  Jacobian  matrix  can  be  formed  by  taking  the  appropriate  partial s: 

r  it  &  v/!  COS  0  rn<  fl  c  in  s'! 


£V/! 


»(5) 


sin  0 


*(!) 

*(!) 


COS  0  sin  0 


i/a 

sin  0 


(nC)1/2  cos  0 
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It  is  now  a  simple  matter  to  show  that,  letting  I.  represent  the  vector  of  the 
ith  column  in  J, 

=  hi]  sij  [hj 

where 


i  *  j 
1  =  J') 


Thus,  the  coordinates  are  truly  curvilinear  and  the  scaling  factors  may  be 
obtained  as: 


h  =  (Hi)*/* 

3 


The  distance  equation  becomes: 

(ds)2  =  (n  +  O(dn)2  +  (n  +  §)(<*$)*  +  n?(de)2  . 

Rewriting  parametrically,  with  £  -  a 

s  (n  +  a)(ax)*  -  na(lf )’] 1/2  dx  • 

There  are  a  variety  of  ways  which  can  be  chosen  to  solve  this  equation,  for 
example 

s  =JT3^(i1  +  a)(^0  +T,a(^)] 1  dx  * 

Let 

n  =  af 
0  =  g 

s  *  aJ  [ir(f  *  1)f'!  *  f9'!]  ,/!  dx 

and  since  the  integrand  is  not  dependent  on  g  explicitly,  then  letting  6  repre¬ 
sent  the  integrand, 


7 


31  October  .1988 


4—  -  c 


[j7(f  +  i)f -  fg'J]  ,/! 
f  V*  =  c!  [-^(f  +  l)f'!  +  fg'«] 


[f>  -  fc!]  g'!  =  |^(f  +  l)f 
,  -  Cf'  ff  +  1  1  l/2 

9  '^LfT^J 


f  -  c2 

f  +  1  =  X2f  -  X2c 

f  -  X2c2  +  1 
X2  -  1  ’ 

r  =  ,  _x  ,  [(,  -  l)2Xc!  -  2X(X*ca  *  1)] 

*  ;  .  ■  [x‘c!  -  ca  -  Xaea  -  l]  . 

(X2  -  l)2  L  J 


(X2  -  l)2 
Continue  with 

9  =  9o  1  f  Jr  xdx 
=  6  ±  |  f ‘  xdx 

0  2J  (X2  -  l)2  X*c2  +  1 

=  e  ±  c(c2  +  1)  f - ^ - 

0  J  (X*  -  l)(X2c2  +  1) 
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now 

1  -  1  f  1  a 

(x  -  l)(ax  +  1)  ”  1  +  a  l_x  -  1  '  ax 

so  that  letting 

a  =  c2 

x  =  X2 


r 


,  =  »<±cj-x.dx(_L-).(_f 

Solving  the  integral: 


c2X2 

c2X2  +  1 


dX 


dX  = 


X* 


dX 


Let 

X  =  tanh  u 
dX  =  sech2  udX 


-  J~ tanh2  udu 

-  f  1  -  sech2  udu 

-  (u  -  tanh  u) 

(X  -  tanh"1  X) 


Let 

cX  =  tan  6 
c*X*  *  tan*  8 

dX  =  ~  sec*  Ode 

C 

obtain: 


t] 

i)  • 


(sech*  u  -  1) 
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sec2  0  -  1  d0 


-[§  Jun‘  0d0 

-  J sec2  9  -  1  dfl 

-  (tan  0  -  0)J  *  "  [c  I  cX  '  tan"1  cX  l] 

g  s  0(  ?  c  [(X  -  tanh'1  X  -  X  -  *  tan-1  cX  J 


0  =  0  ±  £c  tanh'1  X  -  tan-1  ex] 

Recapping  the  first  two  coordinate  variables: 
K  s  a 


... .  rx2c*  +  ii 

’-abrrJ 


In  this  example,  one  constant  of  integration,  c,  appears  in  two  coordinate 
variables,  however,  the  curve  is  still  symmetric  about  an  arbitrary  0fl,  as  one 

would  expect  from  a  surface  of  revolution. 

Continuing  for  4>, 

n-w  dn  :  .  k  d0  * 

0  -  h,  ax  *.  *  h.  as  e» 


°  •  ♦, 

tan  *  =  •+ 

D  •  e^ 

which  resolves  as: 

L  x*  -  i  J 


-2aX(c2  +  1) 

(X*  -  D* 

c(c*  +  1)X2 
(X2  -  l)(X2c2  +  1) 
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h,  =  «n  +  £)■/'  = 


and  obtain: 


X2c2  +  1 V/? 


^  =  tan' 


4 _ 

|_(1  +  c2)1/2 


(x2  -  i y/* 


Finally,  for  distance  transversed: 

As  +  a)(S)  +Tia(ax)]  7  dX 

As  ,  a  f l .  .  (»SL±l  +  i\  .  iL2i5l-U)i 

J  [4(x2c2  +  1)  V  x2  -  l  ;  (x*  -  l)* 

t  X2c2  +  1  (  c2(c2  +  l)2!4  \1  l/2  dx 

X2  -  1  V  (x2  -  l)2  (X2c2  +  l)2  )] 

‘-fcamr* 

Jt&w  "Jin? 


As  =  iC?..-.4'.  1)  [" sinh  2  tanh-1  X  -  2  tanh-1  X 


X  w/approprlate  change  of  limits. 
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PREFACE 

This  is  the  second  addendum  to  AFSPACECOM  Technical  Note,  Surface  Navigation 
and  Geodesy,  A  Parametric  Approach,  1  March  1988.  This  addendum  will  address 
two  classes  of  problems  not  previously  covered. 

First,  a  di cuss ion  of  the  surface  overtake  problem  is  presented.  This  is  a 
classical  naval  problem,  whereby  a  submarine  is  sighted  by  an  intercept  craft; 
the  submarine  submerges  and,  at  an  assumed  constant  velocity,  races  away  at 
some  random  angle  8  (measured  with  respect  to  the  initial  line  of  sight  of  the 
intercept  craft).  Assuming  the  velocity  of  the  intercept  craft  to  be  greater 
than  that  of  the  submarine,  what  is  desired  is  the  curve  the  intercept  craft 
must  take  to  eventually  overtake  and  catch  the  submarine  (this  is  a  well  known 
curve).  The  presented  problem,  solved  on  a  local  planar  surface,  also 
addresses  certain  quantities  of  operations  research  interest,  such  as  mean 
distance  to  intercept,  variance,  mean  time  to  intercept,  etc.  The  overtake 
curve  for  a  sperical  surface  is  then  solved  for  and  presented. 

Second?,  the  problem  of  minimum  distance  navigation  on  composite  surfaces  will 
be  addressed;  specifically,  navigation  on  a  sphere  which  contains  a  spherical 
dent  on  its  surface  (an  approximation  of  the  surface  of  the  martian  moon, 
Phobus) . 


THIS  TECHNICAL  NOTE  IS  THE  WORK  OF  DR  K.  CARLTON- WIPPERN  AND  DOES  NOT 
NECESSARILY  REPRESENT  THE  OFFICIAL  VIEW  OF  THE  U.S.  AIR  FORCE. 
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I.  SURFACE  OVERTAKE  PROBLEM 

This  is  a  special  application  of  the  general  curves  of  pursuit  problem 
(for  example,  reference  Differential  Equations  with  Applications1).  In  this 
case,  for  a  given  surface,  an  interceptee  vehicle  is  sighted  by  the  intercep¬ 
tor.  It  is  assumed  the  interceptee  takes  a  minimum  distance  route  away  from 
where  he  was  sighted,  at  some  random  heading  angle. 

a.  Solution  on  the  plane. 


We  set  up  the  coordinate  system  as  follows,  at  time 


the  inter¬ 


ceptor  is  at  point  A,  the  interceptee  at  the  origin  (point  B),  at  a  distance  d 
from  the  interceptor.  The  interceptee  takes  a  linear  route  away  from  the 
origin  at  velocity  vo,  so  that  at  time  t  -  0  the  interceptee  is  at  a  distance 

d  from  the  origin  but  at  an  unknown  (heading)  location. 


The  equation  of  motion  for  the  interceptee  is  (in  polar  coordinates) 


8  =  fi  (unknown) 

0 

r  =  r0  +  v8t  (where  r#  =  d). 


xRitger  and  Rose,  Different lal  Equations  with  Applications,  p.  94. 
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Assuming  the  interceptor  travels  at  constant  velocity  vx  (vt  >  v0),  the  polar 
equation  of  motion  for  the  inteceptor  becomes 
s2  =  r2  +  r282  (s2  =  v2) 


and  for  eventual  overtake  and  capture,  the  interceptor's  r(t)  component  must 
match  exactly  to  that  of  the  interceptee,  thus 

v2  =  v2  +  (r  +  v  t)2  02 
1  0  '0  0 


2  N  l/2 


X.  K-O 

r  v 

0  0 

1  +  —  t 
r 

0 


=  0 


Letting, 


i.fvl 


1/2 


obtain 


a 


and  solving,  this  becomes: 
a|A0  I 


u  =  u  e 
0 


Substituting  back  for  r 


r  =  re 
0 


a  A8  I 


or  parametrically,  the  solution  curve  for  the  interceptor  is 

r  =  r  +  v  t 
0  0 


(choosing  0#  s  0). 
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b.  Statistical  results  of  the  planar  solution. 

With  the  interceptee  heading  angle  (80)  defined  as  a  random  variable,  it  is  a 

reasonable  question  to  ask  what  are  the  mission  statistics  for  the  average 
mission  (eg,  average  distance  and  time  to  intercept,  and  variance  values). 
Borrowing  heavily  from  Parzen2  we  look  at 


for  the  solution  curve,  then  letting  0  be  a  random  variable 


where  r#  £  r  £  roez1Ia  assuming  an  eventual  overtake  no  later  than  at  0  =  2tt. 
Normalization  for  the  distribution  on  r  becomes 

jira 

1  ■ c  r  » 


so  the  density  function  p(r).is 
r 


P(r)  * 


-  < 


1  .  1 
2ira  r 


r  Sri  re 
0  0 


zito 


•'v 


0  otherwise 

Solving  for  the  average  r 

<r>  =J° rP(r)  dr  =  (eziw  -  1) 
Solving  for  the  variance  of  r;  first 

<r*>  =  f  r'p(r)  dr  =  ^  (.*"*  -  1) 
and  the  variance 

j*  s  <r*>  -  <r>* 


becomes 

•Parzen,  Modern  Probability  Theory  and  Its  Applications,  pp.  308-316. 
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With  the  average  radial  distance,  the  time  to  intercept  becomes: 

<r>  * r, +  ve  *  %  *  r  t<r>  -  r.i 

0 

where  t  is  the  average  intercept  time  (vQ,  rfl  and  <r>  known). 
With  tg,  the  average  total  distance  traveled  per  mission  becomes: 


i  =l2_Y!i\[eiir°  -  a  +  2f«)]  . 

2rm\v0  ) 

If  vx  »  vo  (hence  1/a  »  1  ), 


[Note,  ordinarily  the  density  function  would  be  a  function  of  s,  the  total 

distance  traveled,  not  r;  however,  for  r  -  rQe  r  and  s  are  linearly 
related: 
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so 

ds  =  dr(lt?)‘/'  . 

With  s  and  r  linearly  related,  the  developed  statistics  do  not  change.] 


c.  Solution  on  the  sphere. 


In  this  case  it  is  assumed  the  interceptee  takes  a  great  circle  route  at  some 
unknown  heading  away  from  the  origin  (defined  here  as  in  the  last  example, 
the  point  where  the  interceptor  saw  the  interceptee).  Let's  assume  the  origin 
is  at  some  known  latitude  and  longitude  {La,  Lq}.  We  rotate  the  sphere  so 

that  the  origin  is  at  the  north  pole  of  the  primed  system  (we  do  this  using  a 
topocentric  rotation  matrix3). 


M.  = 


!  sin  La  cos  L# 
-sin  L_ 


sin  <  sin  L„ 
a  o 


.cos  L. 


\  cos  L,  cos  L„  cos  L,  sin  L„ 
\  a  o  a  o 


-cos  la  \ 
0 


sin  La  / 


This  rotation  produces  a  figure  similar  to  the  one  shown  below 


E 


Using  one  more  rotation  about  the  primed  polar  axis 


cos  a 


Mj  =|  sin  a 


-sin  a 
cos  a 
0 


’Mueller,  Bate  and  White,  Fundamentals  of  Astrodynamics,  p.  79. 
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this  aligns  the  interceptor  A  along  the  new  (0  =  0)  axis  as  is  shown  below 

E' 


where  (for  a  sphere  of  radius  a)  the  distance  between  B  and  A  is  d  =  a  •  $fl. 

If  B  is  at  the  origin  at  time  t  =  -  (-■  ^  *  )  then  at  time  t  *  0,  B  will  be 

0 

at  some  random  heading,  a  distance  d  =  a  •  0O  away  from  the  origin: 


In  this  coordinate  system,  a  great  circle  is  represented  by  0  =  60  [consider  a 

vehicle  at  the  real  north  pole,  any  great  circle  route  would  be  due  south, 
along  a  line  of  constant  longitude].  Thus  in  this  coordinate  system,  the 
interceptee's  equation  of  motion  is 

0  =  § 

T  ro  a 

Again  assuming  the  interceptor  travels  at  constant  velocity  vt  (v  >  v0),  the 
equation  of  motion  becomes  ‘  1  “ 
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s2  =  p2  +  p2*2  +  p2  sin2*  6 


(and  since  p  =  a) 


v2  =  a2*2  +  a2  sin2*  0* 

In  analogy  to  the  planar  case,  the  *  components  of  the  two  curves  must  match, 
hence: 


v?  =  <  + 


a2  sin2  ^*0  +  ^tji2 


/v2  -  v2  W2  /  V  v 

(— *r~)  *  s1n(*.  +  rt)° 


d0  =  dt  •  esc 


Letting 


(*.+r 


u  *  ♦.  +  r  * 


i»  4  - 1 


-  *:  \  ,/! 


o  d0  =  esc  u  du 


o0  =  c  -  In 


1  +  cos  u 
sin  u 


ln  -Jln-i?,-  =  tt0  -  c 

1  +  cos  u 


sin  u  _  A  o0 
+  cos  u 


sin  u  _  .  u 
+  cos  u  "  1  » 


and  setting  0  -  0  at  t  =  0  obtain: 


*  ~  2  tan 


-‘[tan  e“9] 


7 


1  May  1989  .  . 


TECHNICAL  NOTE 


or,  parametrically 


e(t) 


- 

v0 

- 

In  tan 

4  +  —  t 
o  a 

-  In  tan 

♦o 

2 

T 

with  (4,  8)  solved  for  in  the  rotated  (or  primed)  system, 
original  (unrotated)  system  can  be  easily  generated,  this 
primed  vector 


the  solution  in  the 
by  creating  the 


f  sin  4  cos  8  ^ 


sin  4  sin  8 


» 


COS  4 


J 


then  finding 


P  =  (MXM0) 


- 1 


Interceptor  latitude 


cos  L. 


cos  L 


a 

sin  L 


cos 

sin 


and  longitude  data  is  available  directly  from 


j 


Mission  statistics  (averages  and  variances  for  a  randomly  headed  interceptee) 
can  be  calculated,  albeit  numerically,  as  before. 

d.  Direct  Intercept  problem. 

This  is  the  standard  curve  of  pursuit  problem.,  whereby  the  interceptor  sees 
the  interceptee,  and  pursues  the  interceptee  along  a  curve  such  that  the 
interceptor's  heading  is  constrained  to  be  simultaneous  with  the  line  of  sight 
to  the  interceptee. 

(1)  Solution  on  the  plane. 
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m  , 

Construct  the  coordinate  system  as  follows: 


y 


_ —d 

f 

c 

-i  i  n 

a  lb 

t  =  0 

where  A  Is  the  interceptee,  initially  at  the  origin,  and  B  is  the  interceptor 
initially  located  a  distance  a  from  the  origin  along  the  x  axis.  Starting  the 
clock  at  time  t  =  0,  then  at  time  t  *  6t  we  assume  the  Interceptee  has  trav¬ 
eled  some  distance  6sfl  along  the  y  axis  at  velocity  v0  and,  similarly,  the 
interceptor  has  traveled  some  distance  6s,,  at  velocity  v(. 

Refering  to  the  figure  below: 


y 


In  time  6t  the  Interceptor  has  moved 

^[-Gor*  ■ 


And  the  Interceptee  has  moved 


«<o 


in  the  same  time. 
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Now,  distance  so,  measured  from  point  B  on  the  curve  y  =  f(x),  is  equal  to 
y  -  x  tan  #  *  y  -  x  and  equating  the  interceptor  to  interceptee  movement 


in  time  fit: 


Letting  a  -  and  assuming  a  <  1  (vQ  <  vx) 

ai(y'  xa?)=*[1+(ax )]  / 


-  X 


day 

dx* 


=  a  1  + 


an 


with  initial  conditions: 


=  0 


> 


y(x  =  a)  =  0 


it  is  relatively  easy  to  show  the  solution  curve  is: 


y 


(2)  Approximate  solution  on  the  sphere.  As  a  first  order 
approximation,  we  will  make  the  following  assumptions: 


(a)  The  Interceptee  Is  located  on  the  equator  and  traveling 
due  east  along  the  great  circle  route  [there  is  no  loss  in  generality  to  this 
assumption,  providing  the  Intercepts*  Is  confined  to  travel  on  a  great  circle 
route.  This  being  assumed  true,  then  a  judicious  choice  of  coordinates  (ie, 
rotating  the  sphere)  will  put  the  Interceptee  on  the  equator  as  desired]. 


(b)  The  Initial  position  of  the  Interceptee  is  analogous  to 
that  of  the  planar  case,  in  that  the  Initial  line  of  sight  of  the  inter¬ 
ceptor  is  perpendicular  to  the  line  of  travel  of  the  Interceptee. 


*R1tger  and  Rose,  Differential  Equations  with  Applications,  p.  512. 
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(c)  The  initial  displacement  of  the  interceptor  from  the  inter- 
ceptee  is  small  so  that  certain  trigonometric  approximations  can  be  made  [eg, 
spherical  position  of  the  Interceptor  {in  #,  0  coordinates}  is  such  that 
sin  f  «  1,  this  being  true  when  the  interceptor  position  is  close  to  the 
equator;  eg,  i  =  90°  ±  20°].  Continuing  as  before,  at  some  time  fit,  the 
placement  of  interceptor  and  interceptee  will  be  as  follows: 


Now,  for  point  B 


fit 


1  +  sin2 


i 


»/* 


d*  •  a 


where  a  is  the  spherical  radius,  and  measuring  s0  from  point  B 


-^  =  0  +  tan"1 

d 


cos  0  sin  * 


The  second  term  is  obtained  from  spherical  trigonometry 


cos  £  -  sin  $  cos  9 
sin  0  =  sin  £  sin  1> 
cos  0  =  cos  £  sin  <t> 

+  sin  £  cos  0  cos  <p 
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tan  0  -  cos  i  tan  f 

and  from  Carton- Wippern*,  where: 


tan  *  -  sin  * 


Equating  both  sides  at  time  fit: 


i-  [l  +  sin2  #(||  y  ‘/!  d*  =  i-  d  [e  +  tan-1  (sin  t  cos  9  f|) 

Now,  using  the  assumption  from  Id(2)(c)  and  by  assuming  that  the  term  cos  0  f 

remains  small  so  as  to  approximate  tan'1  us  u,  the  equation  of  the  pursuit 
curve  becomes: 

viC1  +0ft)  ]  '  =3?[9  +  cos  *(3?). 

Continuing  on  the  right  side  with: 
d0  ,  d20  .  d8  d20 

-rr  +  COS  $  - -  -  Sin  0  -rr  *  cos  $  — - 

a*  d*2  00  d02 

and  letting 


«  -  r  («  <  i)  . 

1 


n  -  dfl 

p  *3? 


a[l  *  p*3‘/*  =  COS  ♦  gj 


a  sec  0  d$  = 


[1  +  p2]1/2 


In  |  sec  #  +  tan  0  |  =  In  |  p  +  (p*  +  l)1/2  |  +  c 


Letting  p  «  0  at  0  =  *0,  then 


"Carl ton- Wippern,  Surface  Navigation  and  Geodesy,  A  Parametric  Approach, 
pp.  38-43. 
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p  +  (p2  +  l)1/2 


where 


v  _  i  .  ir 


tan  (x) 
tan  (xj 


♦/ 


Solving  for  p,  the  differential  equation  for  0  is 

a? =  ?E^tan  *«)""  *  *tan  *)a  ‘  (tan  '  ^tan 

Thus,  the  integral  equation  for  0  as  a  function  of  *,  assuming  convergence, 
becomes : 


(tan  x)“  -  (tan  xo)a  •  (tan  x)"a^j 


d  <t> 


where,  for  the  interceptor,  0  =  0O,  $  =  $0  at  time  t  =  0. 

To  show  that  the  integral  will,  in  fact,  converge  it  is  necessary  to  change 
variables,  +  -  |  -  L,  and  use  L,  the  conventional  latitude  variable. 


Then 


and  the  integral  becomes 


I 


where  as  the  interceptor  approaches  the  interceptee,  L  (or  U)  -*•  0.  In  this 
limit,  it  is  easy  to  show  that  the  above  integral  (in  fact,  the  entire  spheri¬ 
cal  analysis)  converges  to  the  planar  case,  on  the  local  topocentric  plane. 
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II  OVERVIEW  OF  MINIMUM  DISTANCE  NAVIGATION  ON  COMPOSITE  SURFACES 

In  this  section,  we  will  confine  our  attention  to  solving  for  the  curve  of 
minimum  distance  between  two  points  when  the  surface  on  which  the  solution  is 
sought  can  be  described  as  two  simple  surfaces  (eg,  a  sphere  cut  in  half,  so 
that  the  surface  is  a  hemisphere  bounded  on  top  by  a  plane).  We  let  Sx  be  the 

first  simple  surface  with  geodesic  ,  $2  be  the  second  simple  surface  with 
geodesic  $2  and  let  c  =  c(u,  X)  be  the  locus  of  points  which  make  up  the  bound¬ 
ary  between  the  two  surfaces.  If  the  geodesic  between  two  points  traverse 
the  boundary  c  at  say  (ufl,  X0),  then  *x(u0,  Xfl)  =  *2(u0,  Xo)  by  continuity. 
Futhermore,  if  the  boundary  is  smooth  and  continuous, 

(u0»  ^0)  > 

the  tangent  must  also  be  continuous.  While  these  conditions  are  necessary, 
they  are  not  sufficient.  The  problem  is  a  nonlinear  one,  requiring  judicious 
evaluation. 

4  4  4 

Example.  Let  P  ,  P2  and  P3  be  three  points  on  the  equator  of  a  sphere  equally 

•4  *4  *  -4 

spaced  (ie,  Px  at  0°  longitude,  Pg  at  120°  longitude  and  P3  at  240° 

■4  4 

longitude).  The  geodesic  (*x)  between  Px  and  P2  lies  along  the  equator,  from 

■4  -4  -4 

Px,  towards  the  east.  Similarly,  the  geodesic  (*2)  between  P2  and  P3  lies 
along  the  equator,  from  P2  towards  the  east.  Patching  the  two  minima  curves 

4  -4-4-4 

and  *2  together  at  P2  it  Is  observed  that,  at  Pf,  #x(Pt)  =  *2(P2)  and 
T,(Pt)  =  T2(P2), 


■+  3*. 

T.  =  ar 


+  3^, 

(U0»  K)  ~  ^8  ®  3s“ 


•4  -4 

but  the  true  geodesic  from  Px  to  P3  is  not 
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■f  4  *4 

The  geodesic  from  P  to  P3  starts  at  P.x  and  goes  westward,  around  the  back  to 

4 

P  .  The  solutions  for  geodesics  on  more  complicated  surfaces  (such  as  a  com¬ 
posite  surfaces)  can  become  very  involved  and  are  not  necessarily  unique.  In 
the  next  section  we  will  address  a  simple  yet  practical  problem,  solving  for 
the  curve  of  minimum  distance  on  a  sphere  with  a  spherical  dent  feature  on  its 
surface.  Such  a  representative  surface  does  exist  (approximately)  in  nature; 
an  example  would  be  the  martian  moon  Phobus. 

Ill  A  PHOBUS  PROTOTYPE  COMPOSITE  SURFACE  NAVIGATION  PROBLEM 

Consider  the  following  diagram: 


sphere  of  radius  R 


Thus,  the  spherical  dent  in  the  outer  sphere  of  radius  R  is  a  removed  partial 
sphere  of  radius  a.  A  useful  c'oordinate  system  is  one  which  describes  the 
boundary  between  the  two  surfaces  as  a  constant. 
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The  outer  sphere  has  its  center  at  the  origin,  the  boundary  between  the  two 
simple  surfaces  at  4  s  0  .  The  removed  sphere  having  its  center  at  A  and 

corresponding  boundary  at  *  =  $fl.  Letting  the  outer  sphere  be  Sx  with  geode¬ 
sic  and  the  removed  (inner)  sphere  be  S2  with  geodesic  ^2,  five  situations 

can  arise  for  any  two  distinct  points  (PQ  and  Px): 

4  4 

a.  P„  and  P.  may  both  reside  in  S  , 

0  1  *  1  7 

-4  -4 

b.  PQ  and  Px  may  both  reside  in  S2,  or 

4  -4 

c.  P0  and  Px  in  $x  and  $2  respectively  (or  vice  versa). 

-4  -4 

d  &  e.  Po  or  Px  or  both  may  reside  on  the  boundary  between  &  S2. 

* 

Irrespective  of  the  case,  the  boundary  between  Sx  and  S2  not  being  smooth 
requires  us  to  consider  three  curves  or  paths  when  solving  for  the  curve 

4  -4 

of  minimum  distance  between  points  P0  and  Px,  these  paths  are  *  ,  *2*  and  a 
third  path  along  the  boundary  between  the  two  surfaces  (at  0  =  *o,  call  it 

*4  -4 

<f>3)..  Assuming  Po,  Px  are  unit  vectors,  then  in  the  distance  between 
points  is  ASx  =  R  cos-1  (Po  •  Px),  similarly  in  $2,  AS2  =  a  cos"1  (Po  •  Px). 

Example  1. 

•+  -> 

Let  Po  be  on  the  outer  sphere,  Pt  be  in  the  dent.  Find  the  minimum  distance 

-4 

curve.  Let  P($0,  X)  be  on  the  point  on  the  boundary  where  the  minimum 
distance  curve  passes  through.  The  inner  and  outer  spheres  are  assumed  to  .be 

4  -4  A 

longitude  aligned,  that  is,  P(*o,  X)  2  P(*o,  X).  The  idea  is  to  solve  for  X 
that  minimizes  the  total  distance. 

AS(X)  =  R  cos"1  [^P#  •  P(*#,  X)]  +  a  cos’1  [_  ?t  •  P(*Q,  X)] 
and  set 

il*s)  s  0 

3X  u  * 


*which  are  great  circle  curves  on  the  respective  spheres,  Sx  and  S2. 


16 


1  May  1989 


TECHNICAL  NOTE 


fri  general,  a  minima  curve  between  points  may  cross  the  boundary  twice,  at 

P(*,  X)  and  P($,  u).  From  this,  there  are  several  test  cases  which  must  be 
simultaneously  checked,  solved  for  and  compared  to  find  the  true  minima. 

Example  2. 

4  -4 

Let  P#,  Pt  be  on  the  outer  sphere  (S1).  Under  certain  conditions,  the 
minima  curve  may  go  into  S2,  ride  on  the  border  or  remain  totally  in  S  . 

Sample  geodesics 


*  =  1>1 


*  = 


+  A  +  ilfl 


For  each  case  a  test  distance  is  used  AS  =  R  cos-1  (P  •  P  )  provided  all 
along  0  1 
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AS,  =  R  cos‘£p(  •  X)]  +  a  cos- 1  [?(♦„,  X)  •  P(}(,  u)]  + 

R  cos-  »[?(♦,,  u)  •  P,] 
arid 

3AS2  3ASa 
IF  =  ~lu~  =  0  ; 

AS3  =  R  cos-1[p0  •  P(*0,  X)]  +  R  sin  $o  •  [  X  -  u  ]  +  R  cos"1^^,  u)  •  pj 
3AS.  3 AS. 

IF  s  nr =  0  • 

And  so  forth.  Once  P(4,  X),  P(4,  u),  etc.  are  found  which  satisfy  the  minimum 

4  4 

distance  criteria,  it  is  an  easy  matter  to  go  from  P0  to  P  by  way  of  dividing 

up  the  curve  and  using  the  appropriate  parametric  geodesic  relation.5  For 
navigation  on  the  inside  sphere  (S2)  it  is  convenient  to  think  of  the  S2  south 
pole  as  the  S1  north  pole,  ie,  solve  the  S2  portion  as  though  one  was 
navigating  at  or  near  the  south  pole  of  the  excavated  sphere,  but  report  the 

position  latitude  or  ir  -  *  or  northern  hemisphere  latitude.  In  this  way  the 
concept  of  north  pole  remains  consistent,  but  the  curves  must  be  reported 
separately  with  respect  to  each  simple  surface  (inner  or  outer). 

In  summary,  while  not  an  easy  task,  it  is  certainly  possible  (if  not  triv¬ 
ial)  to  examine  on  a  case  by  case  basis,  potential  geodesics  for  the  minimum 
distance  between  points  on  this  simple  composite  surface,  and  then  solving  for 
the  absolute  minima.  Depending  on  the  positions  of  A  (the  excavated  or 
removed  sphere  center)  and  the  ratio  of  a  to  R,  it  may  become  possible  to  make 

*4  *4 

certain  generalizations  (if  Po  and  P2  are  In  Sf  the  geodesic  curve  may  always 
be  ^2,  and  never  include  or  ^  curves;  also,  if  the  removed  spherical  por¬ 
tion  is  less  than  a  hemisphere,  i  ,  the  border  or  rim  geodesic  vanishes  as  a 
potential  solution)  thus,  each  problem  must  be  evaluated  on  its  own  merits. 

£J- 

K.  Carlton-Wippern,  PhO. 

Operations  Research  Analyst 


•Carlton-Wippern,  op.  cit.,  pp.  38-43. 
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PREFACE 

This  is  the  third  addendum  to  AFSPACECOM  Technical  Note,  Surface  Navigation 
and  Geodesy  a  Parametric  Approach,  1  March  1988.  This  addendum  addresses  the 
derivation  of  the  commonly  used  navigation  curve  called  the  rhumb  line  or 
loxodrome. 

The  loxodrome  is  a  well  known  curve,  constructed  so  that  the  angle  between  any 
meridian  and  the  loxodrome  is  constant  throughout  the  curve.1  Although  there 
are  many  ways  to  derive  this  curve,  the  derivation  presented  here  is  based  on 
analysis  provided  in  the  subject  technical  note2  and  on  well  known  results 
from  general  curvilinear  analysis.  One  further  property  of  the  loxodrome 
should  be  mentioned;  when  plotted  on  a  Mercator  projection  the  loxodrome  be¬ 
comes  a  straight  line. 


4The  curve  of  constant  compass  heading. 

2Carlton-Wippern,  Surface  Navigation  and  Geodesy  a  Parametric  Approach,  pp  4, 
33-35. 
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I  INTRODUCTION  AND  GENERAL  DERIVATION 

Consider  u  ,  u,  and  u3  curvilinear  coordinates  with  Jacobian  scaling  factors 
h  ,  h„  and  h,. 

12  3 

Then 

(ds)2  .  (hjduJ2  +  (hadu2)2  +  (h3du3)2  , 

and  if  u  ,  u2  and  u3  can  be  written  as  a  function  of  X 
-►  du3  .  du2  du3  . 

D  =  h3  +  h2  gj-  e2  +  h3  e3  ,*  where  D  is  the  tangent  (or  heading) 

vector. 

Assume  the  surface  under  consideration  can  be  constructed  by  setting  one  of 
the  curvilinear  coordinate  variables  to  a  constant,  ux=  constant  (eg,  for  a 

sphere  p  =  a). 

The  heading  angle,  measured  off  one  of  the  remaining  principle  axes,  can  be 
determine  by 

.  du3 

h3  gr-  h3  du3 

tan  3 . =  r —  •  g-j—  (by  chain  rule).  (1) 

h  au*  n*  au* 

2  3T 

By  assuming  tan  *  a  constant,  equation  (1)  becomes  the  loxodrome  generating 
differential  equation  which  relates  how  uf  is  linked  with  u3.  This  equation 

will  be  used  throughout  this  addendum.  Equation  (1)  can,  In  addition  to 
generating  a  loxodrome  curve,  be  used  to  generate  a  generalMercator  projec¬ 
tion  of  a  surface. 

II  APPLICATION  TO  THE  PLANAR  SURFACE 
With 

and 

•■/ [(S  )••(*)••-(*)■]  • 

obtain  two  loxodrome  forms  on  the  plane,  both  by  setting  z  =  constant.  For 
the  first  case: 

*Car1ton-W1ppern,  op.cit.  p  4. 
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tan  4  =  jjj 
which  solves  as 

y  =  yo  +  tan  <f>  (x  -  xQ)  (a  straight  line). 

For  the  second  case: 


which  solves  as 

r  =  ro  (an  exponential  spiral). 

The  differential  equation  for  the  second  case  may  appear  to  the  reader  to  be  a 
familiar  relation.2 

The  first  case  solution  (a  straight  line)  is  coincidently  the  planar  surface 
geodesic.  A  similar  analogy  exists  for  the  cylinder  (where  helix  is  both 
cylindrical  geodesic  and  loxodrome). 


Ill  APPLICATION  TO  THE  CONIC  SURFACE 
With 


sin' 


f  d0 


T 


i/a 


dX 


(the  surface  generated  by  setting  <t  =  a,  a  constant),  then 


tan  =  P  sin  a  ^  , 

which  solves  as 

cot</>*sina«(8-8„) 
p  =  p  e  o 

r  r  0 

with  *l>  oriented  as  discussed  in  Surface  Navigation  and  Geodesy.3  Distance 
measurement  along  the  conic  loxodrome  is  accomplished  by 

As  =  S  L1  +  p2  si"2  a(a?  )*] l/  dp 

“  f  Q  +  tan2  tfT]  */2  dp 

2Thomas,  Calculus  and  Analytical  Geoaetry,  p.  372,  equation  7a. 
3Carlton-Wippern,  op.cit.  p  19. 
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by  direct  substitution;  then: 

P 

As  -  sec  $  Ap 

Po 

Singularities  occur  when  j  =  *  (n  ±  i)  with  n  =  0,  ±1,  ±2  .  .  .  at  which 
point  p  =  p0  and  c 

0 

As  =  p„  sin  a  •  A0 
r  o 

8o 

Also  of  interest  is  when  *  =  nw  with  n  =  0,  ±1,  ±2  .  .  .  at  which  point  0  =  0o 
and 


p 

As  =  Ap 

Po 

IV  APPLICATION  TO  THE  SPHERICAL  SURFACE 
a.  Spherical  loxodrome. 

With 

15  *  J^'dO  * sin'  *(af)  ]  '  dx 

and  p  *  a,  the  loxodrome  is  generated  by 
-tan  ~  sin  0  ,* 

from  which 

f  *  -cot  ♦  •  <9  • e.> 


*The  minus  sign  introduced  to  reflect  the  proper  orientation,  p  measured  from 
true  north;  therefore,  d#  becomes  -(d*)  so  as  to  produce  a  left-handed  coor¬ 
dinate  system. 
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In  tan  |  =  -cot  •  (8  -  6fl) 


d  .  -cot#*(0-0  ) 

tan  j  -  tan  e  r  v  o' 


For  distance  measurement  along  the  spherical  loxodrome,  by  direct  substitution 


As  =  a J'  [1  +  tan2  $]*/*  d* 


As  -  a  sec  i>  A* 


Again  .two  important  singularities  arise;  +  -  ir  (n  i  i'jwith  n  =  0,  ±1,  ±2 
.  .  here  i  =  a  \  '  J 


As  =  a  sin  A0 


and  *  *  nir  with  n  *  0,  ±1,  ±2  .  .  where  0*0. 


As  =  a  •  A* 


b.  Position  and  velocity  vectors  measured  in  inertial  space. 


As  »  a  sec  *  A* 


4 


TECHNICAL  NOTE J 


1  June  1989 


changes  in  0  can  be  determined,  then  with 

8  =  0o  +  tan  0  •  |jn  tan  -  In  tan  |  J  (0  decreasing) 

changes  in  6  can  be  determined.  Now 
x  V  /  sin  0  cos  [9  +  Q(t-t)] 
y  =  a  sin  0  sin  [0  +  Q(t-t)] 
z  /  \  cos  0 

where  the  term  ft(t-t)  has  been  included  to  correct  for  the  sphere  rotating  at 
an  angular  velocity  Q,  measured  from  some  inertial  body  centered  coordinate 
system. 

For  the  velocity  vector,  let 

/  cos  0  cos  8  0  -  sin  0  sin  9  0 

cos  0  sin  8  0  +  sin  0  cos  0  0 

-  sin  0  0 

0  is  determined  by  As  =  -a  sec  0  A0  (for  0  decreasing)  -s  =  a  sec  00  and 
letting  i  =  v  (a  constant),  then  -  ^  cos  0  -  0;  and  with  -  tan  0  =  sin  0 

then  using  the  chain  rule  -tan  0  0-  sin  0  0,  and  0  =  j  sTn”0 

Now,  by  direct  substitution: 


d 

3t 


y 

=  a 

l*  1 

d 

Ht 


P 

7 

y 

=  V 

\z) 

\ 

sin  0  cos  0  -  cos  0  cos  0  sin  0 


/ 


sin  0  cos  0 

We  add  to  this  vector  the  effect  of  rotation  and  obtain: 

-cos  0  cos  0  cos  5  -  sin  0  sin  8  \  /  -sin  0  sin  8 

y  ]  =  v  |  sin  0  cos  8  -  cos  0  cos  0  sin  6  +  aQ  I  sin  0  cos  8 

sin  0  cos  0  /  \  0 
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where  0  =  0  +  fl(t  -  t) 


-*•  4 


NOTE:  R  •  R  =  0  on  the  static  spherical  surface,  where 


R  -  dR 

’  K  "  3t  • 


V  APPLICATION  TO  THE  OBLATE  SPHERICAL  SURFACE 

a.  Slightly  oblate  spheroidal  loxodrome. 

We  start  by  introducing  the  oblate  spheroidal  curvilinear  coordinate  system: 

x  \  /cosh  u  sin  v  cos  0 

y  =  5  cosh  u  sin  v  sin  0 

z  /  \  sinh  u  cos  v 

This  is  a  curvilinear  system  which  generates  oblate  spheroids  of  eccentricity 
e  by  setting  u  to  a  constant,"  specifically  u  =  sech-1  e.  Since  u  is  set 
constant  we  will  need  only  find  h2  and  h  for 

3 

-</[«(£)■•*!(£)  T*  ' 

and 

hj  =  cosh*  u  cos*  v  +  sinh*  u  sin*  v 
h*  =  cosh*  u  sin*  v  , 
so  with  cosh  u  *  ~ 

As  =  |  J  [d  -  e!  sin*  v)  (  ^  *  sin1  v  )*]  W‘  dx  . 

Now,  recognizing  $  *  ae  (a,  the  semimajor  axis  of  the  ellipsoid)  and  con¬ 
tinuing  as  before: 

(1  -  e 


___  t  _  sin  v  d0 

‘ “ — :* v  yfi  ’ 33 


"Which  is  equivalent  to  setting  curvilinear  (ut)  constant. 
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(1  -  e*  sin*  »)*/»  dv  .  ,cot  4  .  (e  .  e j 


sin  v 


tanh~ 


‘( 


COS  V 


(1  -  e2  sin2 


v)‘/!  ) 


-  e  tanh" 


‘( 


e  cos  v 


(1  -  e2  sin2  v)1 


T7') 


=  -cot  i>  (0  -  80) 
and  with  e2  «  1, 

sTn^  -  f"  sin  v  dv  8  "cot  ♦  *  (fl  -  e0) 


In  tan 


v  .  e 

?  +  r 


COS  V 


s  -  cot  ♦  •  (e  -  eo) 


(for  v  decreasing) 


Thus,  this  is  the  first  order  correction  for  v  [which  by  the  way  is  not  0  in 
spherical  coordinates,  but  is  related  by 

(1  -  e2)1/2  cot  v  =  cot  $]. 

Distance  measured  along  the  loxodrome  becomes: 

As  »  a  f  [1  +  tan*  *]1/2  •  (1  -  e2  sin2  v)1/2  dv 

=  a  sec  i>  f  (1  -  e2  sin2  v)1/2  dv,  and  with  the  integrand  being 
e2 

approximately  1  -  sin2  v 

^  £av  -  ^Av  -  j  sin  2v^J 


As  a  a  sec 


v 

v. 


Solution  singularities  and  resolutions  are  similar  to  those  on  the  sphere. 

b.  Position  and  velocity  vectors  measured  in  Inertial  space. 

With  As  determining  changes  in  v,  which  then  determines  the  changes  in  0,  each 
successive  v,  0  is  found.  The  inertial  position  is  determined  by: 

a  sin  v  cos  [0  +  Q(t  -  x)] 

a  sin  v  sin  [0  +  Q(t  -  x)] 

b  cos  v 

where  a,  b  are  the  semimajor  and  semiminor  axis,  respectively. 

For  the  velocity  vector,,  let 
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a  cos  v  cos  0  v  -  a  sin  v  sin  0  0 


y  =  a  cos  v  sin  8  v  +  a  sin  v  cos  0  0 
\  z  /  \  -  b  sin  v  v 

As  in  the  spherical  case  v  is  determined  by 
ds  =  -a  sec  tp  (1.  -  e2  sin2  v)1/2  dv 

i  =  -a  sec  ip  (1  -  e2  sin2  v)1/2  v  , 
and  letting  as  before  $  *  v,  then: 

v  =  -7  cos  i p  (1  -  e2  sin*  v)"1/2 

Q 


e  =  Y  .  lilLt  , 

a  sin  v  * 

and  by  direct  substitution 

(x  \  /  -cos  v  cos  0  cos  (1  -  e2  sin2  v)'1/2  -  sin  0  sin  ♦  n 

y  =  v  -cos  v  sin  0  cos  tp  (1  -  e2  sin2  v)-1/2  +  cos  0  sin  *  ] 
z  1  \  (1  -  e2)1/2  •  sin  v  cos  tp  (1  -  e2  sin2  v)'1/2  / 

and  finally,  taking  into  account  rotation  about  the  polar  (or  minor)  axis: 
/  *  \  /  -cos  v  cos  9  cos  *  (1  -  e*  sin*  v)"x/2  -  sin  5  sin  *  \ 

y  «  v  -cos  v  sin  0  cos  ♦  (1  -  e*  sin*  v)'1/2  +  cos  S  sin  ^ 

\  t  /  \  (1  -  e2)1/2  sin  v  cos  ♦  (1  -  e2  sin2  v)"1/2  / 


■sin  v  sin  § 


+  afl  sin  v  cos  9 


where  9  =  0  ♦  Q(t  -  x)  . 
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VI  SURFACES  WHERE  THE  LOXODROHE  AND  GEODESIC  CURVES  ARE  COINCIDENT 
Starting  with 

(ds)2  =  h^duj2  +  h2(dua)2  +  h2(du,)2 

and  then  setting  u  =  constant,  thereby  generating  a  surface  spanned  by 

(V 


As  = 


h3  \  /dU3  '2 

i+^ 


h,  du, 


Now,  letting 

u,  =  9, 

du.  . 

3G7  ‘  9 

and  assuming 

3h, 

3g 

ag 

3hJ 

ag 


=  0 
2  =  0 

^  =  0 
ah. 


and 


ag' 


=  0 


/  a(h2  :  h3) 

written  — - —  =  0 

\  a(9  :  9  ) 


the  equation  which  generates  the  geodesic  becomes: 
h 


h  \<  du, 

'•is; I  =7 


/  du,  \  2 

V* 

=  c 


But  by  definition,  the  loxodrome  is  generated  by 
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h,  du3 


tan  i>  =  |p  =  constant. 


Solving  both  equations  simultaneously  by  direct  substitution  yields: 


h  iS!Lt=  c 

"s  sec  * 
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h  =  c  •  esc  i>  , 

3 

and  with  c,  <f>  constant,  this  requires  that  h,  =  constant. 

3 

The  geodesic  loxodrome  curve  relation  between  u2  and  u3  then  becomes: 
u  =  u“  +  c'1  sin  *  tan  *  I  h„  du,  , 

3  3  Z  Z 

where  u“  and  c"1  are  the  two  constants  of  integration  from  the  (geodesic) 

3 

second  order  differential  equation. 


There  are  surfaces  where  the  geodesic  loxodrome  exists  (eg,  on  the  plane, 
the  geodesic  loxodrome  is  a  straight  line;  on  the  cylinder,  the  geodesic 
loxodrome  is  a  helix).  This  mathematical  phenomena  can  occur  when  by  setting 
ux  *  constant,  to  generate  the  surface,  in  turn  implies  that: 

h,  =  constant 

and 


a(g  :  g') 


This  is  not  to  say  that  the  proposed  process  of  generating  the  geodesic 
loxodrome  is  necessarily  unique  [a(h,  :  h,)/3(g  :  g')  =  0  may,  in  fact,  be 

over  restrictive],  but  that  further  analysis  into  the  phenomena  of  the  geode¬ 
sic  loxodrome  is  beyond  the  scope  of  this  addendum. 

On  a  variety  of  surfaces  (eg,  a  sphere)  a  parametric  geodesic  for  a  specific 
case  may  also  be  a  loxodrome.  For  example,  consider  the  geodesic  between 
two  points  on  a  sphere;  one  point  on  the  equator,  the  other  at  the  north  pole. 
The  great  circle  connecting  the  two  points  is  one  of  constant  heading  or 
direction  (due  north).  In  this  case  the  constant  c  *  0,  forcing  i>  =  0  (or 
180°),  and  the  geodesic  loxodrome  equation  reduces  to: 


where  8  is  the  longitudinal  variable. 
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VII  SUMMARY  AND  THE  TORUS  REVISITED 

The  loxodrome,  while  not  as  efficient  a  curve  as  a  geodesic  for  navigation 
purposes,  is  a  useful  navigation  tool  and  is  relatively  easy  to  construct  on  a 
wide  variety  of  surfaces  (provided  there  exists  a  set  of  curvilinear  coordi¬ 
nates  such  that  a  surface  under  investigation  can  be  generated  by  setting  one 
of  the  coordinates  equal  to  some  constant  value). 

As  a  final  comparison  between  the  tractability  of  solving  for  a  geodesic  or  a 
loxodrome  curve,  we  revisit  the  torus.  Choosing  curvilinear  coordinates  such 
that 


Figure  1.  Cross  section  of  a  torus 

In  this  coordinate  system  ht  =  a,  h2  =  av,  h3  =  a(l  +  v  sin  *),  and  setting  v 
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then  solving  for  0  =  8(*), 


0  s  8,  i  a 
0 


/  (IT 


0  u  (1  +  e  sin  #)[(1  +  e  sin  *)* 
However,  solving  for  the  loxodrome: 

.  _  1  +  e  sin  i  d0 
* - ; — 1  3? 


+  e  sin  * 


d$  =  -cot  f  de 


a  ... .  _ ilil! 

Cl  -  Tty/*  Ld  -  €«)i/* 


=  -cot  4>  (0  -  0O) 


NOTE:  The  torroidal  loxodrome  is  not  only  solvable  for  0  lie, 


-  «>)'/*  *  La  -  €*)*/« J  4 


p  »  0;  -  tan  #  •  — 


but  also  invertable  so  as  to  make  i  *  4(0);  given  *0,  0o  ^ie,  calling 
*  2e  .  ,ftan  I/)*  e1 

(i  -  «b)vb  L(i-«,)v,J 

f-  2  tan"1  j(l  -  •  tan  (*0  -  cot  *  (8  -  60))]  -  ej  j 

and  the  measurement  of  distance  along  the  torroidal  loxodrome  becomes: 

As  s  a  S  [€*  +  (1  *  *  sin  #)* (i? )*]  l/*  d* 

As  *  b  J°  jj.  +  •  gi  ^  J  ^  d*  ^remembering  e  =  |  ^ 


As  »  b  sec  +  A# 
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